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1 .  Introduction 


1.1  Objectives  of  the  study 

•  To  produce  a  state  of  the  art  numerical  model  of  a  large  reach  of  the  Missouri  River, 
capable  of  simulating  high  resolution  spatially/temporally  distributed  results  of  water 
depth  and  velocity. 

•  To  investigate  model  behaviour  through  a  sensitivity  analysis. 

•  To  compare  model  predictions  to  multiple  types  of  measured  data,  internal  to  the  model 
domain,  in  order  to  assess  model  performance  and  the  utility  of  the  validation  data. 

•  To  investigate  the  impact  of  varying  levels  of  bathymetric  data  and  mesh  resolutions  on 
the  model  predictions. 


1.2  Background  Information 

The  reach  of  the  Missouri  River  being  used  for  this  modelling  study  is  that  from  Gavins  Point 
Dam,  South  Dakota,  to  the  gauging  station  at  Masked,  Nebraska  (Figure  1.1).  The  reach  covers 
river  miles  81 1  to  776  making  it  35  miles  or  55  km  long.  The  channel  varies  in  width  between 
300m  and  1200m.  The  channel  slope  is  very  low,  dropping  only  about  12  metres  along  its  55 
km  length,  giving  a  gradient  of  0.02%.  The  bed  material  in  this  channel  is  sand  which  is  fairly 
mobile  but  the  channel  banks  have  been  strengthened  or  stabilised  along  much  of  the  reach. 
There  are  several  islands  along  the  reach  several  of  which  are  permanent.  There  is  one  major 
tributary,  the  James  River,  that  joins  the  main  stem  at  river  mile  800  adjacent  to  the  James 
River  Island. 

The  flow  out  of  Gavins  Point  Dam  is  regulated  to  minimise  the  risk  of  flooding  downstream, 
hence  the  reach  being  modelled  is  very  unlikely  to  attain  out  of  bank  conditions.  The  flow  in 
the  James  River  is  however  naturally  variable.  The  model  of  this  system  can  therefore  remain 
as  channel  only.  River  flow  data  is  available  at  several  points  along  the  reach  on  an  hourly 
basis.  Table  1.1  shows  the  data  availability  and  gauge  locations  (also  Figure  1.1). 


Table  1-1  Gauge  station  location  and  data  availability. 


Gauging  Station 

Distance  from  top  of  reach 
(km) 

Data  Available 

Gavins  Point  Dam,  SD. 

0  \  ' 

Flow  rate  ;  '  ^ 

Yankton,  SD. 

Gayville,  SD. 

Masked,  Neb. 

55 

Stage  '  | 

Scotland,  SD. 

53  km  up  James  River 
froiri  confluence 
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Figure  1-la)  Location  and  (b)  detail  of  the  study  reach  of  the  Missouri  River,  USA.  The 
modelled  reach  is  from  Gavins  Point  Dam  to  Maskell  and  includes  the  gauging  stations  at 
Yankton  and  Gayville. 
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The  flow  from  the  James  River  gauge  at  Scotland,  South  Dakota,  is  routed  using  a  simple  one¬ 
dimensional  kinematic  wave  model  to  the  required  location  on  the  main  model  input  boundary 
when  necessary. 

The  hydraulic  model  used  in  this  study  is  TELEMAC-2D.  This  is  a  high  resolution  space/time 
distributed  hydraulic  model  using  finite  element  methodology.  The  model  is  potentially  capable 
of  fulfilling  all  the  objectives  of  this  study.  The  model  is  described  in  some  detail  in  section  2. 


1.3  Key  Features  of  the  study 

•  The  application  of  a  state  of  the  art  two-dimensional  finite  element  code  for  modelling 
large  scale  fluvial  hydraulics  on  the  Missouri  River. 

•  The  use  of  a  high  resolution  model,  in  both  space  and  time,  along  with  wetting  and 
drying  algorithms  for  representing  moving  flow  field  boundaries  allows  dynamic 
inundation  predictions. 

•  For  the  first  time  model  predictions  are  validated  on  this  scale  in  both  time  and  space 
using  multiple  data  sources.  The  data  sources  are  internal  stage  data  at  two  sites, 
supplied  by  the  USAGE  Missouri  River  District  (MRD),  and  satellite  imagery,  supplied  by 
the  Remote  Sensing  and  Geographic  Information  Systems  Center  at  USAGE  Cold 
Regions  Research  and  Engineering  Laboratory  (CRREL). 

•  The  strength  of  model  validation  is  graded  depending  on  the  data  source  and  model 
prediction. 

•  The  influence  of  bathymetric  data  on  the  model  predictions  is  studied  for  its  potential  in 
aiding  future  data  collection  strategy. 
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2.  The  TELEMAC  modelling  system 

The  TELEMAC  system  has  been  developed  by  the  Department  Laboratoire  National 
d'Hydraulique  (LNH)  at  Electricite  de  France,  Direction  des  Etudes  et  Recherches  (EDF-DER). 
LNH  was  formed  in  1946  to  undertake  studies  for  EDF's  hydroelectric  power  projects  and  to 
solve  the  hydraulic  engineering  problems  of  the  Maritime  Ports  and  Waterways  Authority. 
Since  1965  the  role  of  LNH  has  changed  somewhat.  About  70%  of  LNH's  work  is  now  directly 
for  the  benefit  of  EDF,  split  between  hydraulics  inside  power  station  machinery  and  the  study 
of  environmental  problems  in  rivers,  seas  and  the  atmosphere  brought  about  through  the  siting 
of  such  plants.  The  rest  of  LNH's  work  is  on  behalf  of  other  organisations. 

The  TELEMAC  system  therefore  follows  a  line  of  hydraulic  simulation  codes  originating  from 
LNH  for  the  study  of  environmental  problems,  such  as  flooding  resulting  from  a  dam  break  or 
thermal  emission  into  a  river  or  estuary.  TELEMAC  was  developed  from  the  ULYSSE  code,  a 
2D  finite  difference  system  (Ulysses  being  the  father  of  Telemachus  in  Greek  mythology). 
TELEMAC  is  now  in  its  second  four  year  development  period,  each  period  having  a  budget  of 
118  million  francs  (Hardy,  1997).  TELEMAC  is  a  general  purpose  code  that  is  applicable  to 
many  situations  beyond  the  original  remit  of  its  design,  such  as  natural  floods  on  rivers. 

The  TELEMAC  system  is  a  series  of  computer  programs  utilising  finite  element  techniques  for 
simulating  hydraulic  situations.  The  system  includes  pre-  and  post-processing  components  and 
offers  solutions  in  two-  and  three-dimensions.  It  also  contains  facilities  for  sediment  and 
contaminant  transport  (SUBIEF)  and  sand  transport  (TSEF).  The  full  system  is  further  outlined 
in  Figure  2.1.  All  the  components  of  the  system  have  common  file  formats  and  are  written  in 
the  high  level  language  FORTRAN  77.  This  enables  easy  file  transfer  between  components  of 
the  system  and  allows  model  users  to  modify  parts  of  the  code  as  desired  for  specific 
applications. 

The  TELEMAC  system  contains  both  two-  and  three-dimensional  versions.  The  three- 
dimensional  version  (TELEMAC-3D)  has  advantages  in  many  applications  where  vertical 
velocity  variations  are  important  such  as  small  scale  and  oceanic  studies  but  in  large  scale 
fluvial  applications  depth  averaged  calculations  are  adequate.  Indeed  the  equations  used  in 
nearly  all  fluvial  hydraulic  models  derive  from  the  depth  averaged  equations  known  as  the 
Shallow  Water  or  St.  Venant  Equations,  TELEMAC-2D  being  no  exception.  The  third 
dimension,  though  initially  appealing,  would  simply  add  to  the  computational  and  data 
demands.  Descriptions  of  TELEMAC-3D  can  be  found  in  Hervouet  et  al.  (1994),  Hervouet  and 
Janin  (1994)  and  Hervouet  and  Van  Haren  (1996).  TELEMAC-2D  is  however  used  throughout 
the  work  presented  here  and  is  now  discussed  in  more  detail. 
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Limits  of  Domain 


Figure  2-1  The  TELEMAC-2D  modelling  system  (after  Hervouet  and  Lang,  1995). 
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2.1  TELEMAC-2D 

The  main  features  of  TELEMAC-2D  are  listed  by  Hervouet  etal.  (1994)  as: 

•  structured  or  non  structured  meshing, 

•  use  of  Cartesian  or  spherical  co-ordinates, 

•  subcritical  and  supercritical  regimes  (with  hydraulic  jumps), 

•  various  momentum  source  terms  :  bottom  friction,  wind  stress  atmospheric  pressure, 
Coriolis  force, 

•  turbulence  modelling  (k-epsilon  model), 

•  equation  on  temperature  or  a  substance  concentration, 

•  treatment  of  tidal  flats, 

•  many  types  of  boundary  conditions  including  free  slip  condition  and  incident  waves. 

Many  of  these  are  of  minimal  importance  for  fluvial  applications  but  are  necessary  for  the 
diverse  nature  of  cases  that  TELEMAC-2D  can  be  applied  to.  A  wide  range  of  test  cases  are 
illustrated  in  the  TELEMAC-2D  (version  3.0)  validation  document  (Cooper,  1996)  prepared  to 
substantiate  the  explicit  claims  made  about  the  applicability  and  accuracy  of  the  computer 
code.  The  test  cases  illustrated  range  from  the  Western  European  Coast  and  the  Mersey  estuary 
UK,  to  flows  around  bridge  piers  and  over  a  weir.  On  a  scale  more  relevant  to  this  study  dam 
breaks,  flow  at  a  river  confluence  and  the  simulation  of  a  flood  event  on  the  River  Culm,  UK 
are  also  shown.  A  new  version  of  TELEMAC-2D  is  released  annually.  All  of  the  simulation  in 
this  report  have  been  carried  out  using  version  3.0,  released  in  1995. 

2.1.1  Shallow  water  equations 

TELEMAC-2D  solves  the  shallow  water  equations  (SWE),  the  depth  averaged  version  of  the 
fully  three  dimensional  Navier  Stokes  equations  of  fluid  flow.  The  SWE  require  that; 

•  the  flow  is  incompressible 

•  the  water  column  is  well  mixed,  so  that  there  are  no  significant  density  variations  in  the 
vertical 

•  vertical  accelerations  are  negligible  (hydrostatic  pressure  approximation) 

•  the  effective  lateral  stresses  can  be  represented  by  an  isotropic  (usually  constant)  eddy 
viscosity 

•  bed  stresses  can  be  modelled  using  a  linear  or  quadratic  friction  law. 

All  of  these  conditions  are  commonly  met  in  rivers,  estuaries  and  seas  making  the  choice  of 
these  equations  over  the  full  three  dimensional  Navier  Stokes  equations  valid  for  the 
applications  of  the  model  (Cooper,  1996). 

The  choice  of  formulation  of  the  SWE  used  in  TELEMAC-2D  is  not  obvious.  A  conservative 
form  would  seem  better  but  divisions  by  the  water  depth  are  needed  to  produce  the  velocity 
field,  hence  causing  problems  in  dry  areas  of  the  model  domain.  Hence  a  non-conservative 
form  is  preferred.  Moreover  numerical  stability  analysis  also  favours  the  non-conservative 
version  of  the  equations.  Using  finite  element  methods  mass  conservation  can  be  ensured  with 
non-conservative  equations.  Two  versions  of  these  non-conservative  equations  have  been 
developed,  the  celerity-velocity  version  and  the  depth-velocity  version.  The  depth-velocity 
version  are  marginally  favoured  as  they  are  easier  to  employ  mass-conservation  techniques  to 
(Hervouet  and  Janin,  1994).  These  equations  are  shown  below: 

2L. 

—  +  u.grad(h)+  h.div(u)  =  0 
5t 
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9u  ,  dh  dZf 

—  +  u.  grad(u)  +  g  -  div{v.  grad(u))  =  Sx  -  g — 
dt  dx  ox 


2.2 


dv  dh  .  „  ctZi 

—  +  u.grad(u)  +  g  —  -  div(v.grad(v))  =  Sy  -  g— 
5y  ay 


2.3 


where: 

h  :  water  depth 

u,  V  :  velocity  components 

g  :  gravity 

Zf  :  bottom  elevation. 

Sx,  Sy  :Source/sink  terms  (bottom  friction,  wind,  etc.) 

V  :  eddy  viscosity 

h,  u  &  V  are  the  unknown  variables. 

Equation  2.1  being  the  continuity  equation,  2.2  and  2.3  the  force-momentum  equations.  Full 
derivations  can  be  found  in  Norton  etal.  (1973). 


2.1.2  Solving  the  equations 

The  TELEMAC  system  uses  finite  element  methodology  to  solve  the  shallow  water  equations 
to  produce  values  of  water  depth  (h)  and  two  velocity  components  (u  and  v)  at  all  points  in  the 
model  domain.  To  achieve  this  the  domain  must  be  discretized  into  a  grid  or  mesh,  usually 
made  up  of  linear  triangles  (unstructured  grid)  for  flexibility.  The  mesh  is  created  outside 
TELEMAC  using  specialised  mesh  generation  software  such  as  I-DEAS  or  BALMAT.  The 
mesh  is  then  incorporated  with  the  topographic  data  to  produce  a  geometry  file  that  is  then  used 
in  the  TELEMAC-2D  simulation.  The  mesh  consists  of  a  series  of  node  points,  the  vertices  of 
the  triangles  with  a  elevation  (z)  value,  where  the  equations  are  solved  and  the  triangles 
themselves  are  the  elements.  The  geometry  notation  used  is  shown  in  Figure  2.2. 

All  three  dependant  variables  (h,u  and  v)  are  defined  at  each  point  in  the  domain  in  terms  of 
linear  interpolation  functions  associated  with  the  value  attained  at  each  node.  Linear 
interpolation  functions  between  nodes  can  be  a  limitation  where  the  bed  gradients  are  varying 
rapidly  compared  to  the  distance  between  nodes.  Mesh  refinement  can  however  be  used  to 
overcome  this  problem  in  most  cases.  The  time  discretization  is  semi-implicit  so  a  system  of 
3N  simultaneous  algebraic  equations  is  obtained,  where  N  is  the  number  of  nodes  in  the  domain 
(Cooper,  1996).  Furthermore  TELEMAC-2D  employs  element  by  element  (EBE)  techniques  to 
minimise  and  simplify  computation.  More  detailed  descriptions  of  finite  element  methodology 
can  be  found  in  Norton  et  al.  (1973)  and  Finder  and  Gray  (1977). 
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Several  solution  algorithms  are  available  in  TELEMAC-2D  including  (Cooper,  1996): 

•  a  fractional  step  method  using  characteristics  for  the  advection  terms  and  a  Galerkin 
method  for  the  propagation  and  diffusion  terms  in  the  equations 

•  several  variants  of  the  Streamline  Upwind  Petrov  Galerkin  (SUPG)  method  (Brooks  and 
Hughes,  1982) 

•  a  hybrid  scheme  that  combines  the  characteristics  and  SUPG  approaches. 

Further  details  of  these  options  are  available  in  Hervouet  and  Van  Haren  (1995). 

In  the  simulations  presented  later  in  this  report  a  fractional  step  method  (Marchuk,  1975)  is 
used  where  the  advection  terms  are  solved  initially  followed  by  a  second  step  where  the 
propagation,  diffusion  and  source  terms  are  solved.  For  the  advection  step  the  method  of 
characteristics  is  used  for  the  momentum  equation  and  a  SUPG  method  for  advection  of  h  in 
the  continuity  equation,  ensuring  mass  conservation.  The  use  of  the  semi-implicit  SUPG 
method  for  the  continuity  equation  achieves  unconditional  stability  (Hervouet  and  Janin,  1994). 
The  second  step  of  the  calculation,  propagation,  is  then  solved  using  a  conjugate  gradient  type 
method. 

The  different  solution  techniques  do  however  produce  very  similar  results  when  used  with  on 
the  same  simulations  with  the  same  mesh.  For  example  Bates  et  al.  (1996)  show  the  SUPG  and 
hybrid  methods  used  on  a  flood  event  on  the  River  Culm,  UK,  showing  very  similar  outflow 
hydrographs  and  inundation  extent. 

The  Courant  number  is  used  as  a  measure  of  the  quality  of  the  numerical  solution.  It  derives 
from  explicit  numerical  methods  and  is  calculated  thus: 


C,=u  — 

Ax 


where  (Bates  et  al.  1996): 
:  Courant  number 


t 

X 

u 


:  time  step 
:  mesh  size 
:  flow  velocity. 
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Explicit  methods  become  unstable  when  C>\.  Although  TELEMAC-2D  is  implicit  and 
theoretically  not  subject  to  such  constraints  the  Courant  number  remains  a  good  measure  of  the 
quality  of  the  numerical  solution  (Hervouet  and  Janin,  1994;  Bates  et  a/.,  1996).  Hervouet  and 
Janin  (1994)  suggest  that  in  some  cases  simulations  can  be  performed  with  values  of  the 
Courant  number  of  up  to  50  but  this  is  not  advisable. 

2.1.3  Boundary  conditions 

There  are  two  main  types  of  boundary  condition  that  can  be  used,  solid  boundaries  and  liquid 
boundaries.  The  boundaries  described  here  are  round  the  sides  of  the  model  domain  and 
through  the  bed  of  the  model.  All  boundary  conditions  are  assigned  on  a  node  by  node  basis. 
Solid  boundaries  are  no  flux  (impermeable)  and  incorporate  a  friction  factor  (Hervouet  and  Van 
Haren,  1995).  Liquid  boundaries  allow  a  flux  across  them.  They  are  more  difficult  to  deal  with 
as  they  suppose  the  existence  of  a  fluid  area  that  is  not  part  of  the  calculation  domain  but  can 
however  influence  it.  This  influence  is  described  through  the  boundary  condition.  There  are 
four  types  of  liquid  boundary,  entry  and  exit  with  supercritical  flow  (Froude  number  >1)  and 
entry  and  exit  with  subcritical  flow  (Froude  number<l).  Incident  waves  and  prescribed 
flowrates  can  be  incorporated  through  these  4  boundary  types.  Hervouet  and  Van  Haren  (1995) 
describe  these  boundary  conditions  more  fully. 

2.1.4  Physical  parameter  options 

Physical  representation  of  parameters  is  essential  in  models  such  as  TELEMAC-2D  in  order  to 
apply  them  to  different  applications.  TELEMAC-2D  includes  numerous  physical  parameters. 
The  most  important  are  discussed  in  this  section. 

2. 1. 4. 1  Bottom  friction 

There  are  six  options  in  TELEMAC-2D  for  representing  bed  friction,  these  are: 

•  No  friction 

•  Linear  friction 

•  Chezy's  law 

•  Strickler's  law 

•  Manning's  law 

•  Nikuradse's  law 

The  applied  friction  coefficient  is  converted  to  force  terms  in  the  x  and  y  directions  at  each 
computational  node  which  is  then  fed  into  the  Shallow  Water  Equations  (in  the  source  terms  Sx 
and  Sy  -  see  section  2.1.1.,  equations  2.2  and  2.3).  In  reality  the  bed  friction  force  is  a  quadratic 
function  of  velocity  so  the  no  friction  and  linear  law  are  rarely  used  in  practical  applications 
(Hervouet  and  Van  Haren,  1996).  The  Chezy,  Strickler  and  Manning  laws  are  all  closely  related 
and  utilise  the  quadratic  function  mentioned  and  are  all  described  more  fully  by  Hervouet  and 
Van  Haren  (1995).  Nikuradse's  law  calculates  a  Chezy  coefficient  from  the  water  depth  (h)  and 
grain  size  of  the  bed  material  which  is  then  converted  into  a  force  term.  Which  of  these  laws 
that  is  used  is  not  very  important  as  they  are  very  closely  related  and  friction  coefficients  can 
easily  be  converted  between  them.  The  choice  rests  with  the  model  user. 
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2. 1. 4. 2  Turbulence  representation 

Turbulence  is  one  of  the  remaining  major  unresolved  problems  in  physics.  Its  importance  in 
hydraulics  is  well  known  but  it's  representation  in  hydraulic  models  is  a  different  problem.  The 
concern  in  hydraulic  modelling  is  how  turbulence  in  flow  affects  the  mean  structure  of  the 
flow.  Turbulence  modelling  is  a  scientific  discipline  in  itself  and  the  reader  is  referred  to 
Younis  (1996)  for  a  general  discussion  of  turbulence  modelling  and  to  Hervouet  and  Van  Haren 
(1996)  for  a  fuller  exposition  of  turbulence  representation  in  TELEMAC-2D  than  is  possible 
here.  In  TELEMAC-2D  version  3.0  there  are  two  options  for  turbulence  representation: 

•  constant  viscosity 

•  k-epsilon  model. 

The  first  and  simplest  is  a  constant  viscosity  coefficient,  zero-equation  turbulence  model.  The 
term  VELOCITY  DIFFUSIVITY  is  used  to  set  the  molecular  viscosity,  turbulent  viscosity  and 
dispersion  throughout  the  model  domain.  In  turbulent  flows  the  turbulent  viscosity  is  dominant 
and  can  virtually  be  equated  to  the  velocity  diffusivity  value.  A  constant  value  of  turbulent 
viscosity  over  the  model  domain  is  often  used  but  is  an  oversimplification. 

The  second  is  the  k-epsilon  model  where  the  turbulent  viscosity  is  expressed  as  a  function  of 
the  turbulent  kinetic  energy  (k)  and  its  dissipation  rate  (s).  This  two-equation  turbulence  model 
is  overcomplicated  for  large  scale  applications  and  is  computationally  intensive. 

A  third  turbulence  model,  the  single-equation  Elder's  model,  has  just  been  introduced  into  the 
newest  version  of  TELEMAC-2D  (version  3.1).  This  model  separates  longitudinal  and 
transverse  dispersion  terms.  It  offers  improvements  on  the  previous  two  methods  in  fluvial 
applications  where  such  separation  is  important.  Unfortunately  this  method  was  not  available 
for  the  work  in  this  report. 

2.1.43  Other  physical  parameters 

Several  other  physical  parameters  can  also  be  specified  in  TELEMAC-2D.  These  include  wind 
stress  on  the  water  surface,  atmospheric  pressure,  water  temperature  and  water  density.  The 
Coriolis  force  can  also  be  applied  when  modelling  large  areas. 

2.1.5  Wetting  and  drying  zones 

Applications  of  TELEMAC-2D  to  river  and  estuaries  generally  involve  wetting  and  drying 
areas  of  the  model  domain.  For  example  the  tidal  simulation  in  the  Mersey  estuary  and  the 
flood  event  on  the  River  Culm,  both  described  by  Cooper  (1996),  involve  this  type  of 
behaviour.  The  ability  of  the  model  to  deal  with  this  sort  of  behaviour  is  therefore  of  vital 
importance  but  is  problematic.  The  behaviour  in  the  dry  zones,  where  divisions  by  the  water 
depth  (h)  in  the  calculations,  can  cause  spurious  terms  to  appear  as  h  tends  towards  zero. 

Two  solutions  to  this  problem  are  available  in  TELEMAC-2D,  both  described  in  more  detail  by 
Hervouet  and  Van  Haren  (1995), 

•  solving  the  equations  everywhere  and  coping  with  the  spurious  terms. 

•  removing  the  dry  zones  from  the  computational  domain. 

The  first  is  the  simplest  but  corrections  must  be  applied  in  the  dry  zone.  In  dry  areas  the  water 
surface  gradient  becomes  that  of  the  bottom  topography  but  this  cannot  be  allowed  to  act  as  a 
driving  force  in  tbe  momentum  equation. 
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The  second  method  removes  the  dry  zone  from  the  computational  domain  and  is  often  called 
the  "moving  boundary  technique".  In  TELEMAC-2D  this  is  achieved  efficiently,  avoiding  the 
need  to  redefine  the  mesh  at  each  time  step,  by  keeping  the  elements  in  the  mesh  but  cancelling 
their  existence  through  the  use  of  an  array  set  to  0  for  dry  elements  and  1  for  all  others. 

Partially  dry  elements  are  another  important  area,  especially  where  the  elements  are  large  such 
as  in  the  Missouri  River  model.  These  are  coped  with  in  TELEMAC-2D  by  a  sophisticated 
method  that  allows  the  water  depth  to  go  to  zero  within  an  element.  This  compares  favourably 
with  the  usual  techniques  of  keeping  the  elements  fully  wet  or  excluding  partially  dry  elements 
(Figure  2.3). 


2.2  Post-processing  of  results 

After  the  computations  have  been  completed  the  numerical  results  must  be  converted  to  a  more 
user  friendly  visual  form.  This  is  done  using  the  graphics  visualisation  software  RUBENS,  also 
written  by  EDF-DER.  RUBENS  allows  quick,  easy  and  flexible  representation  of  the  results  in 
many  different  forms  such  as: 

•  mesh  plots 

•  vector  fields 

•  contour  lines 

•  coloured  surfaces 

•  space  profiles 

•  time  profiles 

RUBENS  also  allows  the  superposition  of  measurements  and  graphics  and  also  the 
manipulation  and  processing  of  the  visualised  data.  Most  of  the  plots  of  model  results  in  this 
report  are  created  using  the  RUBENS  software.  For  more  information  on  RUBENS  the  reader 
is  referred  to  Piro  (1993). 


2.3  Summary  of  the  TELEMAC-2D  system 

TELEMAC-2D  is  a  high  resolution  space/time  distributed  hydraulic  model  that  solves  the 
Shallow  Water  Equations  for  fluid  flow  using  finite  methodology.  The  model  can  be  used  in  a 
wide  variety  of  scenarios  including  those  involving  wetting  and  drying  fronts  within  the  model 
domain.  Bed  friction  and  turbulence  are  represented  in  the  model  through  the  use  of  physically 
based  parameters. 

TELEMAC-2D  can  therefore  be  seen  to  be  well  specified  for  the  type  of  fluvial  application  that 
is  of  interest  in  this  study.  The  code  has  been  shown  to  work,  through  the  report  of  Cooper 
(1996),  in  a  wide  range  of  cases.  Success  in  any  individual  situation  is  however  dependant  on 
the  data  provided,  physical  and  numerical  parameters  used.  How  these  important  factors  have 
been  determined  in  the  Missouri  River  model  case  is  the  subject  of  the  next  section. 
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Figure  2-3  Partially  dry  elements:  (a)  the  real  situation,  (b)  element  remaining  fully  wet 
introducing  spurious  momentum  terms,  (c)  element  exclusion  method  and  (d) 
TELEMAC-2D  method  keeping  the  element  partially  dry  (after  Price,  1997). 
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3.  The  Missouri  River  model 

This  section  looks  at  how  TELEMAC-2D  has  been  applied  to  the  Missouri  River  between 
Gavins  Point  Dam,  South  Dakota,  and  Masked,  Nebraska.  The  theoretical  aspects  of  the  model 
were  examined  in  section  2  but  the  practical  details  of  applying  the  model  to  this  specific 
scenario  reach  is  vital  to  gain  an  understanding  of  the  model  and  its  capabilities. 


3. 1  Producing  the  finite  element  mesh 

The  finite  element  mesh  for  the  model  must  be  created  by  first  defining  the  boundary  of  the 
model  domain  and  secondly  discretizing  this  area  into  elements. 

The  boundary  of  the  2D  model  was  defined  by  digitising  round  the  edge  of  the  river  as 
represented  on  United  States  Department  of  the  Interior  7.5  minute  quadrangle  series  maps. 
The  use  of  wetting  and  drying  algorithms  in  the  model  enable  the  flow  field  boundary  to  be 
calculated  within  this  outer  boundary  thus  allowing  this  boundary  to  be  fairly  loosely  defined. 
The  digitised  boundary  was  then  converted  to  Universal  Transverse  Mercator  (UTM)  metre  co¬ 
ordinates  for  compatibility  with  the  metric  scale  used  by  TELEMAC-2D.  The  boundary  of  the 
model  (with  TELEMAC-2D’s  metric  scale),  which  includes  three  permanent  islands,  is  shown 
in  Figure  3.1. 


A  mesh  of  triangles  or  quadrilaterals  must  be  made  within  the  model  domain  to  facilitate  the 
finite  element  solution  technique.  The  finite  element  meshes  used  with  the  Missouri  model 
were  generated  inside  the  prescribed  boundary  using  the  mesh  generation  package  BALMAT. 
The  meshes  created  were  made  of  near  equilateral  triangles  in  order  to  increase  the  accuracy 
and  minimise  mass  conservation  errors.  Two  meshes  have  been  created  of  this  reach.  The  two 
are  shown  in  Figure  3.2  and  their  attributes  compared  in  Table  3.1. 
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Mesh  1,  the  lower  resolution  mesh,  is  used  for  the  majority  of  the  simulations  in  this  report.  It 
should  be  assumed  this  has  been  used  when  looking  at  results  unless  otherwise  stated.  Mesh  2 
is  only  used  in  the  sensitivity  analysis  (section  4)  and  section  6  where  its  performance  is 
compared  to  that  of  mesh  1 . 


Table  3-1  Comparison  of  the  attributes  of  the  two  meshes  used  to  model  the  Missouri 
River  in  this  report. 


Attribute 

Mesh  1 

Mesh  2 

Number  of  Nodes 

.  .5969  - 

■  9213  . 

Number  of  Elements 

■  10437 

16567 

Average  Element  Size  (m^) 

<  7938  \ 
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3.2  Bathymetric  Data 

The  bathymetric  data  in  a  2D  hydraulic  model  is  one  of  the  most  important  factors  in  attaining 
high  quality  simulations.  The  quantity  and  quality  of  the  bathymetric  data  is  therefore  of  utmost 
significance. 

Along  this  reach  of  the  Missouri  River  the  bathymetry  was  obtained  by  the  MRD  in  1995  by 
echo  sounder  surveys  of  channel  cross  sections.  In  total  343  cross  sections  were  available  along 
the  55  km  reach  creating  an  average  spacing  of  just  over  150m.  Each  cross  section  was  made  up 
of  between  30  and  200  data  points,  depending  on  section  length.  The  typical  spacing  between 
data  points  along  a  section  was  6m.  All  elevation  data  was  converted  to  metres.  Bathymetric 
data  was  missing  on  the  cross  sections  where  the  bed  was  above  the  water  surface,  such  as  on 
sand  banks  and  mud  flats.  In  these  regions  the  topographic  elevation  was  estimated  from  United 
States  Department  of  the  Interior  7.5  minute  quadrangle  maps,  except  on  the  permanent  island 
in  the  model  where  no  such  data  is  required. 

The  bathymetric  data  must  now  be  applied  to  the  finite  element  mesh.  This  is  the  process  of 
interpolating  between  topographic  data  points  and  assigning  elevation  (z)  values  to  the  nodal 
points  in  the  mesh.  This  produces  a  geometry  file  for  the  model  simulation.  The  interpolation  is 
usually  carried  out  using  STBTEL,  the  TELEMAC  sub-program.  This  uses  a  quad-directional 
interpolation  routine,  taking  the  nearest  data  point  in  the  four  quadrants  around  each  node, 
weighting  for  distance  from  the  node  and  combining  to  give  the  nodal  elevation.  An  alternative 
interpolation  routine  has  been  written  at  Bristol  for  use  with  cross  sectional  data  of  river 
channels.  This  is  a  linear  interpolation  down  the  line  of  greatest  depth  between  sections  and 
improves  the  definition  of  this  line  in  such  cases  as  the  Missouri  River.  A  portion  of  the  models 
topography  is  shown  in  Figure  3.3. 
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Figure  3-3  An  example  portion  of  the  models  bathymetric  representation.  This  section  is 
immeadiateiy  upstream  of  Yankton. 
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3.3  Boundary  Condition  Specification 

The  boundary  conditions  of  the  model  are  very  important  to  the  simulations  run.  There  are 
several  boundaries  to  each  model  that  must  be  set  with  care  prior  to  the  model  run.  These 
boundaries  are  divided  into  two  types,  liquid  boundaries  and  solid  boundaries. 

With  the  Missouri  model  the  liquid  boundaries,  allowing  flow  across  them,  are  at  the  top  and 
bottom  ends  of  the  reach  and  on  the  James  River  tributary  inflow.  Model  inflows  are  prescribed 
as  flowrates,  at  Gavins  Point  Dam  and  from  the  James  River,  and  outflows  as  water  surface 
elevations,  at  Maskell.  This  produces  a  well  posed  problem  for  fluvial  flows  according  to  the 
theory  of  characteristics  (Hervouet  and  Van  Haren,  1995).  All  measured  flowrates  were 
converted  to  m^/s  (cumecs)  and  stage  values  to  metres. 

Solid  boundaries  allow  no  flux  across  them.  They  are  found  down  the  sides  of  the  reach  and 
through  the  bed.  The  side  boundaries  in  the  Missouri  model  are  set  as  slip  boundaries,  allowing 
a  velocity  along  them.  This  is  justified  by  the  size  of  the  elements  making  flow  overly 
restrained  in  areas  of  channel  constriction  when  the  more  realistic  no-slip  boundaries  are 
imposed.  The  flexible  boundary  of  the  flow  field  removes  this  problem  for  large  portions  of  the 
reach  where  such  an  argument  is  irrelevant. 


3.4  Physical  Parameter  Specification 

The  two  most  important  physical  parameters  in  fluvial  hydraulic  models  are  generally  agreed  to 
be  bed  friction  and  turbulence  (Baird  and  Anderson,  1990;  Bates  et  al,  1992).  The  theory  of  the 
two  in  TELEMAC-2D  has  been  discussed  in  section  2.1.4.  They  are  the  only  two  considered  in 
any  detail  here. 

Given  the  channel  only  nature  of  the  model  and  lack  of  additional  information  the  bed  friction 
and  turbulence  parameters  were  defined  as  constant  throughout  the  entire  reach.  These  are 
obviously  simplifications  but  should  be  adequate  as  a  first  approximation.  The  bed  friction  is 
always  defined  using  Manning's  law  enabling  the  standard  Manning's  'n'  measure  of  flow 
resistance  to  be  used.  The  turbulence  was  defined  using  the  zero-equation  velocity  diffusivity 
representation. 


3.5  Creating  the  initial  conditions 

Before  any  dynamic  simulation  can  take  place  the  model  must  be  run  to  a  position  closely 
resembling  those  required  at  the  start  of  the  simulation,  the  initial  conditions. 


3. 6  The  model  run 

The  model  run  is  the  end  product  of  the  modelling  process.  This  run  simulates  the  flow  event 
that  is  required  creating  the  model  predictions.  In  order  to  carry  out  this  run  all  the  previously 
discussed  factors  must  be  in  place.  A  geometry  file  containing  the  topography  and  finite 
element  mesh,  a  boundary  condition  file,  parameter  values  and  the  initial  conditions  file  must 
all  be  accessed.  A  time  step  for  the  model  run  must  also  be  set.  Given  the  implicit  nature  of  the 
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numerical  schemes  in  TELEMAC-2D  the  size  of  the  time  step  should  be  relatively 
unimportant.  It  does  however  influence  the  Courant  number  and  therefore,  despite  probably  not 
causing  instability,  a  large  time  step  could  undermine  the  quality  of  the  simulation.  With  the 
Missouri  River  model  the  time  step  used  is  usually  4  seconds  creating  Courant  numbers  less 
than  unity. 


3. 7  Summary  of  model  application 

The  application  of  TELEMAC-2D  to  the  Missouri  River  between  Gavins  Point  Dam  and 
Masked  has  been  described  in  detail.  The  data  sets  are  all  of  high  quality  and  the  model  set  up 
is  fairly  simple  but  justifiable.  The  model  is  now  in  a  position  to  be  applied. 
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4.  Sensitivity  Analysis 

The  first  stage  in  model  application  is  the  sensitivity  analysis.  Sensitivity  analysis  is  a  widely 
used  technique  in  hydrological  and  hydraulic  modelling  to  determine  the  impact  of  changing 
parameter  values  and/or  input  stresses  on  the  model  predictions.  Virtually  every  modelling 
study  involves  a  sensitivity  analysis  at  some  stage.  It  is  a  technique  that  is  potentially  useful  in 
model  formulation,  model  calibration  and  model  verification  (McCuen,  1973). 


4. 1  Background  to  the  sensitivity  analysis 

Howes  and  Anderson  (1 988)  suggest  that  a  sensitivity  analysis  can  be  used  to: 

•  Demonstrate  that  in  response  to  representative  variation  of  model  input  and  parameter 
values,  theoretically  realistic  model  behaviour  is  experienced. 

•  Illustrate  the  model  to  be  sufficiently  sensitive  to  represent  actual  variation  in  the  prototype 
system. 

•  Identify  those  model  parameters  or  inputs  to  which  the  model  is  most  sensitive. 

•  Assess  model  behaviour  without  recourse  to  comparisons  to  field  data. 

Perhaps  more  importantly  in  the  case  of  this  model  is  that  sensitivity  analysis  is  very  useful  for 
deciding  which  parameters  to  adjust  during  calibration.  Adjusting  the  most  sensitive  parameters 
will  produce  a  greater  improvements  in  the  predictions  for  small  changes  in  parameter  values, 
which  means  that  if  the  initial  values  were  approximately  physically  realistic  then  the 
calibrated  ones  should  be  as  well.  Conversely  if  the  less  important  parameters  are  left  fixed  at 
incorrect  values  then  the  resulting  error  on  model  predictions  is  likely  to  be  small  (Troutman, 
1985). 

Some  of  the  theoretical  considerations  of  sensitivity  analysis  are  now  discussed.  Sensitivity  is 
defined  as  the  rate  of  change  in  one  factor  with  respect  to  change  in  another  factor. 
Mathematically  it  can  be  derived  from  a  Taylor  series  expansion  of  model  behaviour  and  after 
discarding  the  non-linear  terms  the  linearized  sensitivity  (S)  can  be  expressed  in  simple  terms 
as  (McCuen,  1 973) 

S  =  3F,/3Fi  4.1 

where  Fo  the  model  output  prediction  and  F|  is  the  model  input  parameter.  In  many  instances 
this  linear  approximation  of  sensitivity  is  adequate  to  describe  model  behaviour  but  whether  it 
is  appropriate  for  the  Missouri  model  is  unknown  at  present. 

There  are  two  methods  for  determining  the  value  of  the  model  sensitivity  which  are  (Lane  et 
al,  1 994). 

•  The  direct  differentiation  of  the  governing  equations.  This  has  been  demonstrated  with 
simple  models,  for  example  Beven  (1 979)  uses  this  method  to  assess  the  sensitivity  of  the 
Penman-Monteith  evapotranspiration  equations  and  LaVenue  et  al.  (1 989)  use  it  to  help 
estimate  travel  time  uncertainties  in  ground  water  flow.  However  with  complex,  distributed 
models  this  methodology  has  not  been  developed  sufficiently  and  cannot  possibly  assess 
the  entire  complex  response  of  a  model  such  as  TELEMAC-2D. 

•  Factor  perturbation  approaches  can  evaluate  the  sensitivity  of  the  model  by  incrementing 
parameter  values  and  assessing  the  model  response.  This  is  the  most  commonly  used  form 
of  sensitivity  analysis. 
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Given  the  problems  of  direct  differentiation  in  assessing  the  sensitivity  in  complex  models  the 
parameter  perturbation  approach  is  the  only  one  considered  in  this  rest  of  this  section. 

Assuming  the  overall  physically  based  validity  of  TELEMAC-2D  to  this  type  of  application  the 
sensitivity  analysis  of  the  Missouri  River  model  has  one  major  function.  This  to  set  up  the  basis 
for  the  following  calibration  procedures  by  determining: 

•  the  relative  sensitivity  of  the  different  physical  parameters, 

•  the  effect  of  mesh  resolution  on  this. 

Once  these  questions  have  been  resolved  the  application  and  calibration  of  the  model  can 
proceed. 

The  runs  carried  out  in  this  sensitivity  analysis  are  outlined  in  Table  4.1.  The  data  from  these 
simulations  is  enough  to  fulfil  the  two  aims  above.  As  discussed  in  the  previous  section  only 
two  physical  parameters  in  the  model  were  varied,  bed  friction  and  turbulence  representation 
(velocity  diffusivity),  both  of  which  were  evenly  distributed  across  the  domain. 


Table  4-1  The  simulations  carried  out  for  the  sensitivity  analysis.  The  figures  in  bold  are 
the  parameter  values  held  constant  whilst  the  other  is  varied. 


Parameter  Varied 

Mesh  1 

Mesh  2 

Bed  Friction  (n) 

0.01,0.015,  0.02,  0.025, 

0.03,  0.035,  0.04 

Turbulence  (Velocity 
Diffusivity  -  m^/s) 

1,2,4  • 

The  range  of  bed  friction  was  chosen  to  represent  the  very  broad  range  for  this  channel  type, 
sand-bed  with  a  shallow  gradient  (0.02%  in  this  case)  (Table  4.2).  This  full  range  is  used  here 
as  no  further  data  is  available  to  make  a  further  judgement  at  this  stage  and  uncertainties  in  the 
overall  model  set-up  create  the  possibility  that  the  values  needed  to  fit  the  results  to  the 
observed  are  not  those  that  would  be  selected  from  the  field.  The  velocity  diffusivity  values  are 
harder  to  physically  justify  given  their  vastly  oversimplified  application  to  the  flow.  The  values 
chosen  are  similar  to  those  used  in  other  applications  that  have  performed  well,  producing 
realistic  output  and  maintaining  model  stability. 


Table  4-2  Values  of  Manning’s  n  for  various  channel  types  (after  Bathurst  1988). 


Channel  Type 

Channel  Slope  (%) 

Manning  ‘n’  range 

Sand-bed 

Gravel-bed 

Boulder-bed 

0.03  -  0.2  '  - 

Steep  pool/fall 

.  ■  >5, 

-  0.1-5^' 

The  parameter  space  was  sampled  using  simple  single  parameter  perturbation  techniques.  The 
flow  record  used  for  these  sensitivity  tests  is  that  from  around  the  6th  June  1994  when  the  river 
is  essentially  at  a  steady  state.  The  results  are  therefore  all  steady  state  values  of  the  variables. 
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4.2  Results  of  the  sensitivity  analysis 

The  results  of  the  parameter  changes  are  looked  at  on  the  inundated  area,  water  depth  at 
Gayville  and  Yankton,  velocity  at  Gayville  and  the  spatial  distribution  of  differences  of  water 
depth  and  velocity  between  two  runs  over  the  whole  domain.  At  this  stage  no  comparison  to 
observed  data  is  made. 

Firstly  the  relative  sensitivity  of  the  two  parameters  on  both  meshes  is  looked  at.  Figure  4.1 
shows  the  impact  of  parameter  variation  on  the  total  percentage  inundation  extent  in  the 
domain.  It  can  clearly  be  seen  from  this  that  the  bed  friction  has  a  far  greater  impact  than  the 
velocity  diffusivity  over  the  ranges  tested. 
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Figure  4.2  shows  the  impact  of  both  bed  friction  and  velocity  diffusivity  on  the  water  surface 
elevation  at  both  gauge  stations,  Yankton  and  Gayville,  for  both  mesh  resolutions.  The  results 
in  all  cases  show  that  the  bed  friction  has  a  far  greater  influence  on  the  water  surface  elevation 
than  the  velocity  diffusivity. 
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Figure  4-2  Sensitivity  analysis  plots  for  water  surface  against  bed  friction  at  (a)  Yankton 
and  (b)  Gayville  and  against  velocity  diffusivity  at  (c)  Yankton  and  (d)  Gayville. 


The  impact  of  parameter  variation  on  velocity  at  a  point,  in  this  case  at  Gayville,  is  shown  in 
Figure  4.3.  Bed  friction  is  again  dominant  but  the  response  is  not  linear. 
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Figure  4-3  Change  in  velocity  at  Gayville  as  (a)  bed  friction  and  (b)  velocity  diffusivity 
varies  for  both  high  and  low  resolution  meshes. 


From  the  above  three  diagrams  it  is  clear  that  the  bed  friction  is  clearly  the  most  influential 
parameter  having  a  far  greater  impact  than  velocity  diffusivity  on  all  the  model  results.  The 
mesh  resolution  has  a  consistent  impact  on  the  model  predictions.  The  greater  the  resolution  of 
the  mesh  the  lower  the  prediction  for  an  equivalent  parameter  set  for  all  model  results.  This  is  a 
slightly  unusual  result  and  shows  that  the  model  structure  does  have  an  important  impact  on  the 
model  predictions. 

Having  determined  that  the  bed  friction  is  by  far  the  dominant  parameter  in  the  Missouri  River 
model  the  further  analysis  is  carried  varying  only  this  parameter.  Using  the  argument  of 
Troutman  (1985),  the  velocity  diffusivity  can  only  introduce  a  small  error  into  the  results  and 
varying  it  produces  only  small  changes,  hence  it  can  reasonably  be  overlooked.  The  difference 
in  velocity  and  water  depth  were  calculated  at  all  node  points  between  the  highest  friction  run 
(n  =  0.04)  and  the  lowest  friction  run  (n  =  0.01)  on  the  Mesh  1.  The  results  are  then  plotted  up 
as  a  spatial  plot  of  the  domain. 
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The  water  depth  difference  plot  (Figure  4.4a)  clearly  shows  that  along  the  thalweg  the  model 
behaviour  is  expected  in  that  the  water  depth  is  greater  with  higher  friction  (positive  difference 
values).  On  areas  with  very  little  or  no  water,  such  as  sand  banks  in  the  chaimel  or  out  of 
chaimel  regions  in  the  model,  the  behaviour  is  very  much  less  marked  with  only  very  small 
changes  being  present,  as  would  be  expected.  There  are  however  areas  in  the  model  domain 
where  the  behaviour  looks  to  be  more  complex  and  perhaps  unexpected. 


The  velocity  difference  plot  (Figure  4.4b)  looks  to  more  complex  in  structure  than  the  water 
depth  one.  Along  the  thalweg  the  velocity  is  generally  lower  with  the  higher  friction  (negative 
difference  values)  but  in  the  shallower  regions  the  opposite  is  true.  This  is  due  to  the  water 
depth  increase  in  these  shallow  areas  facilitating  an  increase  in  velocity  despite  the  increased 
bed  friction.  There  are  however  many  areas  in  the  model  where  the  behaviour  is  unexplained. 
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Further  analysis  of  these  plots  has  been  carried  out  to  find  out  more  about  the  complex  model 
behaviour  at  all  the  nodes.  The  nodal  values  of  velocity  difference,  water  depth  difference,  and 
the  water  depth  for  the  low  friction  run  were  extracted  from  TELEMAC-2D  results  files.  The 
velocity  difference  was  then  plotted  against  the  two  water  depth  variables  mentioned.  The 
results  are  of  great  interest  (Figure  4.5)  and  show  that  unexpected  model  behaviour  is  much 
more  frequent  than  anticipated.  Figure  4.5a  shows  the  velocity  difference  plotted  against  the 
water  depth  in  the  low  friction  run  and  in  this  case,  from  the  above  discussion,  it  would  be 
expected  that  at  all  nodes  but  those  of  very  shallow  depth  the  velocity  difference  should  be 
negative.  This  is  however  clearly  not  the  whole  story.  Although  the  majority  of  velocity 
differences  are  negative  there  are,  across  the  whole  depth  spectrum,  nodes  that  show  the 
opposite  behaviour.  At  very  low  (near  zero)  depths  the  majority  of  nodes  do  show  an  increase 
in  velocity  but  this  is  not  complete.  The  magnitude  of  the  velocity  differences  seems  to  be 
greater  at  smallish  depths  (around  Im)  possibly  caused  by  bed  friction  having  a  greater 
influence  on  flow  hydraulics  at  such  depths  than  at  greater  depths.  Although  with  the  limited 
number  of  points  at  greater  depths  this  could  just  be  a  random  effect. 
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Figure  4-5  Sensitivity  analysis  plots  showing  results  for  all  node  points  for  differences 
between  high  and  low  friction  runs  for  (a)  velocity  difference  against  water  depth  in  the 
low  friction  run  and  (b)  velocity  difference  against  water  depth  difference. 


Figure  4.5b  shows  the  velocity  difference  against  the  water  depth  difference.  The  first  point  to 
note  is  that  there  are  points  that  unexpectedly  show  a  decrease  in  water  depth  as  friction 
increases.  All  these  points  also  show  a  decrease  in  velocity.  Between  water  depth  changes  of 
approximately  0  to  0.5m  there  is  predominantly  an  increase  in  velocity  as  friction  increase, 
probably  corresponding  to  the  nodes  with  shallow  water  depths  initially.  At  greater  water  depth 
differences  the  majority  of  the  behaviour  is  as  expected,  decreasing  velocity  as  friction 
increases  but  there  is  a  significant  proportion  of  nodes  revealing  the  opposite,  though  not  to  the 
same  magnitude. 

These  results  have  shown  the  complexity  of  the  model  response  to  very  simple  changes  in 
parameter  values.  Some  of  the  unusual  results  can  be  explained  by  physical  reasoning,  others 
are  perhaps  artefacts  generated  by  the  model  specification. 


4.3  Conclusions  of  the  sensitivity  analysis 

The  sensitivity  analysis  has  shown  that  the  bed  friction  is  by  far  the  most  influential  parameter 
on  all  outputs  in  this  model  application.  The  strategy  for  calibrating  the  model  must  therefore 
be  based  on  varying  this  parameter.  At  present  the  parameters  are  distributed  evenly  over  the 
entire  model  domain  but  although  this  is  perhaps  simplistic  there  is  no  case,  or  data,  to  change 
them  at  the  moment.  The  mesh  resolution  has  a  consistent  influence  on  the  model  predictions 
showing  the  model  structure  is  an  important  factor  in  the  determining  the  model  predictions 
and  thus  diminishes  the  physical  basis  of  the  model.  The  sensitivity  analysis  also  highlighted 
some  strange  and  unexpected  non-linear  res{5onses  in  some  areas  of  the  model  domain. 
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5.  Comparison  of  model  predictions  to  observed  data 

The  TELEMAC-2D  modelling  system  has  now  been  shown  to  be  a  high  quality  hydraulic 
model  capable  of  simulating  flows  in  a  variety  of  situations.  The  setting  up  of  the  model  and  its 
subsequent  sensitivity  analysis  have  shown  the  Missouri  River  application  to  be  a  highly 
specified  and  physically  based  hydraulic  model.  The  final  tests  for  the  model  are  the 
comparison  of  model  predictions  to  observed  data.  The  observed  data  for  the  Missouri  model  is 
of  very  high  quality  and  quantity.  The  available  data  sets  are  of  special  value  because  they 
allow  both  detailed  point  comparisons  using  the  stage  data  at  two  sites,  Yankton  and  Gayville, 
and  spatially  distributed  comparisons  using  the  satellite  imagery.  Using  this  data  the 
performance  of  the  Missouri  River  model  can  be  comprehensively  analysed. 

The  comparison  of  model  predictions  to  field  data  is  a  common  task  in  all  hydrologic 
modelling,  as  such  a  broad  literature  exists  on  this  subject.  It  is  still  however  a  very  difficult 
area.  The  latest  theories  on  this  are  that  models  can  never  be  completely  validated  for  two 
reasons.  Firstly  hydrological  models  are  essentially  complex  hypotheses  and  can  therefore 
never  be  proved  correct  only  supported  or  disproved  (Konikow  and  Bredehoeft,  1992;  Beven, 
1993).  Secondly  models  can  be  considered  as  a  series  of  equations  with  more  variables 
influencing  the  behaviour  than  equations  to  account  for  them.  The  models  cannot  therefore 
account  for  all  possible  behaviour  and  hence  cannot  be  validated.  Instead  of  complete 
validation  it  is  now  the  task  to  look  for  partial  or  achievable  validation  of  the  model  (Samper  et 
al,  1990).  As  such  the  strength  of  the  partial  validation  can  be  gauged  depending  on  which 
aspect  of  model  performance  the  data  validates.  As  was  indicated  in  the  introduction  the 
determination  of  the  strength  of  validation  of  the  2-D  hydraulic  model  with  the  novel  data 
sources  is  one  of  the  objectives  of  this  project. 

The  observed  data  used  must  be  comparable  to  a  model  prediction.  These  can  therefore  be 
water  depth  (free  surface  elevation)  and  its  two  dimensional  extension  inundation  extent,  flow 
velocity  and  the  outflow  hydrograph.  The  strength  of  these  can  be  ranked  according  to  three 
criteria,  where  >  indicates  greater  strength  of  validation: 

•  internal  >  external  -  whether  the  validation  data  is  internal  or  external  to  the  modelled 
domain. 

•  2D  >  1 D  -  whether  the  validation  data  is  1  -  or  2-dimensional.  For  a  2-D  model  the  2-D  data 
is  the  more  powerful  as  it  cannot  be  simulated  by  a  simpler  model  and  is  more  difficult  to 
match. 

•  process  >  numerical  -  process  validation  is  matching  the  processes  occurring  in  the  field 
with  those  in  the  model,  this  is  usually  considered  more  powerful  than  simple  numerical 
matching  of  results  that  could  have  been  reached  by  incorrect  processes  (Fawcett  et  al, 
1995). 

The  data  available  for  this  study  is  stage  data  at  two  sites  within  the  model  domain  which  can 
be  considered  internal  but  1-D  and  numerical.  The  satellite  images  are  also  internal  but  2-D  and 
numerical  (except  for  example  of  process  validation  that  shall  be  shown  later).  The  satellite 
imagery  is  therefore  the  more  powerful  validation  data  type.  How  well  the  observed  and 
predicted  match  shall  be  shown  in  the  following  sections  and  at  the  end  this  idea  of  strength  of 
validation  shall  be  revisited  in  the  light  of  the  results.  Firstly  however  the  methods  of  model 
application  used  in  this  study  and  elsewhere  are  reviewed. 

Usually  when  hydrological  models  are  applied  to  actual  scenarios  a  two  stage  process  of  model 
calibration  and  validation  is  carried  out.  Calibration  is  the  process  by  which  parameter  values 
are  varied  within  reasonable  ranges  until  the  differences  between  observed  and  computed 
values  are  minimised  (Konikow  and  Bredehoeft,  1992).  Theoretically  speaking,  physically 
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based  model  should  not  need  calibrating.  The  estimation  and  application  of  physically  realistic 
parameter  values  for  the  model  should  enable  the  model  to  perform  well  without  any  further 
manipulation.  However  the  numerous  simplifications  involved  in  modelling  and  problems  in 
parameter  estimation  mean  that  such  a  notion  is  unrealistic.  Calibration  is  therefore  performed 
on  virtually  all  physically  based  hydrologic  models.  Following  the  calibration  phase  the  model 
must  then  be  validated.  In  most  cases  the  calibrated  parameter  values  are  used  again  on  a 
different  portion  of  the  flow  record  to  see  if  they  enable  a  good  match  between  observed  and 
predicted  variables.  If  they  do  then  the  model  and  parameters  are  considered  validated  and  can 
be  used  with  confidence  for  prediction  of  future  events.  If  not  then  calibration  and  possibly 
some  stage  of  model  formulation  must  be  repeated  until  the  model  is  validated.  It  should  be 
noted  at  this  point  that  validation  is  used  here  to  mean  that  a  reasonable  comparison  is  obtained 
between  observed  and  predicted  data  rather  than  that  the  model  has  been  shown  to  be  a  accurate 
representation  of  the  system. 

With  this  application  of  TELEMAC-2D  the  highly  structured  methodology  for  calibration  and 
validation  outlined  above  is  not  entirely  appropriate  for  several  reasons: 

•  this  study  is  for  research  of  model  behaviour, 

•  the  number  and  variety  of  data  sources  make  its  application  difficult, 

•  there  is  no  need  to  produce  a  model  capable  of  predicting  future  events  at  this  stage. 

The  two  processes  are  therefore  not  used  as  distinct  entities  in  the  following  investigation. 
Instead  a  hybrid  of  the  two  forms  is  used  enabling  all  benefits  of  the  calibration-validation 
procedure  to  be  accrued  but  also  much  more.  The  model  shall  be  run  using  a  wide  range  of  bed 
friction  values  over  several  sections  of  the  flow  record.  This  shall  enable  the  calibrated 
parameter  values  to  be  found  for; 

•  each  comparison  data  set  (2  internal  stage  gauges  and  possibly  some  areas  on  a  satellite 
image) 

•  each  portion  of  the  flow  record  utilised. 

Cross  comparison  of  parameter  values  and  the  assessment  of  errors  when  the  model  is 
calibrated  onto  one  data  source  allow  an  assessment  of  the  model/parameter  performance  more 
complete  and  powerful  than  with  the  traditional  calibration-validation  procedure. 

The  sections  of  the  flow  record  to  be  utilised  have  very  different  flow  conditions  and  occur 
over  a  decade  meaning  that  bed  morphology  changes  could  have  occurred.  Three  of  the 
sections  of  record  coincide  with  satellite  images  of  the  reach.  The  flow  records  around  the  dates 
of  the  images  showed  very  little  (<5%)  variation  in  the  data  set  enabling  the  runs  to  remain  as 
steady  state  simulations,  ideal  for  comparing  to  the  snapshot  images.  One  other  period  of 
steady  state  flow  shall  be  utilised,  being  shortly  after  an  image  but  with  different  flowrates  it 
should  enable  a  comparison  between  the  two  sections  of  the  flow  record  without  any  impact  of 
bed  morphology  change.  Finally  one  dynamic  section  of  the  flow  record  shall  be  simulated. 
Table  5.1  shows  all  the  flow  records  to  be  simulated  and  the  available  data  for  comparisons. 


Table  5-1  The  flow  records  to  be  simulated. 


Event 

Number 

Dates 

Type 

Inflows  -  Gavins  Point  ; 
James  River  (m^/s) 

Satellite  Image 

Event  1 

6th  June  1994  : 

Steady  State 

,995  :;t5- 

LANDSAT-TM 

Event  2 

25th  June  1994 

Steadj^Stllte 

Event  3 

26th  June  1984 

Steady  State 

LANDSAT-TM 

Event  4 

8  th  February  1991 

Steady  State 

Event  5 

Dynamic  . 

fQ5  -  765:4S  liil 

None  available 
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The  following  analysis  is  split  between  the  different  portions  of  the  flow  record.  Each  is  dealt 
with  separately,  comparing  the  available  data  and  finding  the  calibrated  parameter  values 
before  all  the  data  is  brought  together  again  for  an  overall  comparison.  The  analysis  of  the  first 
event  contains  details  of  the  background  and  methods  used  in  that  and  all  the  other  cases. 


5.1  Event  1  (6th  June  1994) 

The  time  period  being  simulated  for  this  comparison  is  around  6th  June  1994  synchronous  to 
the  satellite  image  of  the  area  that  will  be  used  in  the  following  section.  The  flow  record  around 
this  date  showed  very  little  (<2%)  variation  in  the  data  set.  The  simulation  was  therefore  taken 
as  a  simple  steady  state  computation,  simplifying  the  analysis  of  results.  All  the  following 
comparisons  are  therefore  done  on  a  single  point  in  time  basis,  e.g.  a  single  water  level  for  eaeh 
point  in  the  domain  for  eaeh  parameter  set. 

5.1.1  Internal  Stage  Data 

Two  internal  gauges  for  water  level  are  available  for  the  reach  of  the  Missouri  River  described 
in  section  1.  The  gauges  are  at  Yankton,  7.5  km  down  the  reach,  and  Gayville,  21.5  km  down 
the  reach  (Figure  1.1).  Hourly  stage  records  are  available  for  these  gauges. 

The  sensitivity  of  the  model  prediction  of  the  water  depth  at  these  two  locations  has  been 
shown  to  be  dependant  predominantly  on  the  bed  friction  parameter.  This  is  the  only  parameter 
varied  in  the  following  analysis.  Other  parameters  would  have  only  minimal  effect  but  increase 
complexity  greatly.  The  sensitivity  of  water  depth  to  friction  is  virtually  linear  at  both  internal 
gauge  sites  (Figure  5.1)  and  in  the  expected  direction,  i.e.  water  depth  increasing  as  friction 
increases. 


This  finding  allows  the  use  of  McCuen's  (1973)  linearized  sensitivity  equation  seen  in  section  4 
but  a  simple  regression  equation  is  more  useful  for  applying  the  results  in  this  analysis  and  is 
therefore  used  here.  Equations  were  produced  for  the  model  behaviour  at  both  sites,  Yankton 
and  Gayville,  relating  the  error  in  water  level  prediction  to  the  bed  friction  (Manning's  'n') 
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parameter.  The  equations,  calculated  using  the  statistical  computing  package  MINITAB,  can  be 
written  thus: 

Yankton  Error  =  -0.410  +  34.72  n 
Gayville  Error  =  -0.976  +  3 1 .93  n 

Alternatively  they  can  relate  the  water  depth  to  the  bed  friction  thus: 

Yankton  Water  Depth  =  4.87  +  34.72  n 
Gayville  Water  Depth  =  2.91  +  31.93  n 

The  confidence  in  these  regression  lines  is  huge  with  tiny  residuals  (<0.01 5m)  in  all  cases. 
These  equations  can  therefore  be  used  as  a  very  powerful  tool  for  estimating  model  predictions 
for  'n'  values  at  which  the  model  has  not  been  run,  although  it  should  not  be  extrapolated 
beyond  the  range  of 'n'  values  used  to  make  the  relationship  (i.e.  0.01  to  0.04).  At  present  this  is 
only  shown  at  two  points  but  could  perhaps  be  applied  at  every  node. 

It  can  easily  be  calculated  from  the  above  equations  that  the  'n'  value  to  eliminate  the  error  in 
water  surface  elevation  at  Yankton  is  0.01 18  whereas  at  Gayville  it  is  0.0306.  The  power  of  the 
above  regression  equations  can  now  be  utilised  to  calculate  the  error  in  the  other  observation 
whilst  one  is  correct.  Whilst  the  stage  at  Yankton  is  predicted  correctly  at  Gayville  the 
prediction  is  0.60m  too  low.  Reversing  this,  whilst  Gayville  is  predicted  correctly  the 
prediction  at  Yankton  is  0.65m  too  high.  In  percentage  terms  these  errors  are  18.26%  and 
10.96%  of  the  expected  water  depth  at  the  two  points  respectively. 

Taking  the  value  central  to  the  above  two  estimates  of  Manning's  'n'  should  enable  the 
approximate  calculation  of  the  joint  minimum  errors  in  the  predictions.  This  value  of  'n'  is 
0.0212  and  produces  an  error  at  Yankton  of  0.33m  and  at  Gayville  of  -0.30m.  Slight  adjustment 
of 'n'  could  equalise  these  errors  but  would  be  irrelevant.  These  errors  can  be  expressed  again  as 
percentages  of  the  water  depth  such  that  the  error  at  Yankton  is  5.89%  and  at  Gayville  is 
8.37%.  Figure  5.2  shows  a  downstream  section  of  the  comparison  between  observed  and 
predicted  water  surface  values. 
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Interestingly  the  value  of  rate  of  change  of  the  water  surface  (error  term  or  water  depth)  is 
different  at  the  two  gauging  stations.  Thus  a  fixed  change  in  the  bed  friction  produces  a 
different  change  in  the  water  surface  elevation  at  the  two  stations.  Investigation  into  how  this 
type  of  behaviour  occurs  over  the  whole  reach  would  be  useful  in  determining  more  about 
model  behaviour  or  what  causes  this  phenomenon. 

The  comparison  between  observed  an  predicted  stage  values  has  shown  that  the  model  is 
capable  of  predicting  this  variable  to  about  0.3m  at  two  sites  simultaneously.  This  is  a  good 
result  given  the  simplicity  of  the  model  set  up.  Given  further  work  on  the  parameterization  this 
could  no  doubt  be  reduced  significantly. 

5.1.2  Satellite  Imagery 

Many  types  of  satellite  imagery  have  been  used  widely  in  studies  of  fluvial  and  coastal 
hydraulics.  The  wide  coverage,  ease  of  use  and  wealth  of  information  contained  in  the  images 
makes  them  ideal  for  many  uses  in  this  field.  Redfern  and  Williams  (1996)  review  the  available 
sensors  and  their  applications  in  (primarily  coastal)  hydraulics.  Previously  a  lot  of  work  has 
been  carried  out  on  the  remote  sensing  of  flood  extent  because  of  the  immense  hazard  and 
mitigation  costs  on  floodplain  developments.  Rango  and  Salomonson  (1974)  show  that  the 
areal  extent  of  flooding  can  be  mapped  using  near-infrared  sensors  on  ERTS  1  (Earth 
Resources  Technology  Satellite  -  later  renamed  Landsat).  More  recently  Imhoff  et  al.  (1987) 
use  SAR  (Synthetic  Aperture  Radar)  and  Landsat  MSS  (Multispectral  Scanner)  for  mapping 
flood  extent  and  damage  in  Bangladesh.  SAR  has  been  shown  to  very  accurate  for  flood 
boundary  delineation  by  Biggin  and  Blythe  (1996)  on  the  River  Thames  in  the  UK.  Satellite 
remote  sensing  of  the  1993  floods  on  the  Mississippi  and  Missouri  Rivers  in  the  USA  with  both 
SAR  (Brackenridge  et  al,  1994)  and  Landsat  TM  (Thematic  Mapper)  have  further  shown  the 
potential  of  such  methodology. 

Satellite  data  has  also  been  used  to  measure  flood  stages  as  illustrated  by  Koblinsky  et  al 
(1993)  using  the  U.S.  Navy’s  Geosat  radar  altimeter  on  the  Amazon  basin  and  Brackenridge  et 
al  (1994)  using  SAR  data  and  topographic  maps  on  the  Mississippi.  Both  methods  have  quite 
large  errors  associated  with  them  at  present. 

Plumes  of  sediment  rich  water  can  also  be  identified  as  they  have  a  higher  reflectance  than 
clear  water  in  the  visible  region  of  the  spectrum  (Lillesand  and  Kiefer  1987).  Brackenridge  et 
al  (1994)  use  this  factor  to  highlight  levee  breaches  along  the  Mississippi  River  during  the 
1993  floods,  the  breaches  acting  as  a  sediment  sources,  and  explain  sediment  deposition 
patterns.  Other  forms  of  pollution  and  thermal  emissions  can  also  be  traced  (Lillesand  and 
Kiefer,  1987). 

The  availability  of  satellite  imagery  for  the  modelled  reach  of  the  Missouri  River  synchronous 
to  the  flow  records  has  enabled  several  fundamental  research  objectives  to  be  addressed.  Data 
from  satellite  images  is  ideal  for  use  with  2-D  fluvial  hydraulic  models  being  of  a  scale 
commensurate  to  the  model  resolution  and  being  widely  spatially  distributed.  In  this  study  such 
data  has  the  potential  to  be  used  for; 

•  validation  of  inundation  extent  predictions, 

•  identify  regions  of  error, 

•  calibration  data, 

•  process  validation. 

Remotely  sensed  data  has  the  potential  to  be  used  more  widely  in  hydraulic  modelling  for 
topographic  and  physical  parameterization  (Bates  et  al,  1997)  but  these  are  ongoing  research 
themes. 
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The  images  available  for  this  study  are  Landsat  TM  for  events  1  and  3  and  SPOT  (Systeme 
Probatoire  pour  I'Observation  de  la  Terre)  for  event  4.  The  resolution  of  Landsat  TM  is  30m 
and  SPOT  is  10m,  both  higher  than  the  resolution  of  the  model  so  they  are  both  more  than 
adequate  to  supply  the  quality  of  data  needed  to  achieve  the  above  objectives. 

Six  sites  have  been  chosen  to  be  looked  at  along  the  modelled  reach.  They  have  been  chosen 
specifically  because  of  their  complex  bathymetry  and  topography,  involving  mud  flats, 
permanent  and  temporary  islands.  The  location  of  the  six  sites  are  shown  on  Figure  5.3.  Side  by 
side  comparisons  of  the  observed  and  predicted  inundation  are  to  be  made  using  different 
friction  values  at  all  six  sites  followed  by  an  overlay  of  the  two  showing  the  geo-referencing 
that  can  be  done  to  the  results. 


X  :zzzz  zzzzz  zzzzz  ^zzzz 


Figure  5-3  Locations  of  the  six  sites  to  be  used  in  the  comparison  of  model  predictions 
with  satellite  imagery. 


The  side  by  side  comparisons  of  the  model  predictions  against  the  Landsat  TM  image,  using 
bed  friction  values  of  0.01,  0.02,  0.03  and  0.04,  are  shown  in  Figure  5.4a-f  for  the  six  sites. 
These  figures  plot  the  satellite  image  against  the  predicted  flow  field  boundary  as  this  is  the 
model  result  that  can  be  directly  compared  to  the  image.  From  these  it  can  be  seen  that  the 
model  predictions  are  generally  good,  producing  a  close  spatial  match  against  the  observed 
data.  There  are  some  obvious  areas  where  the  bathymetry/topography  in  the  model  is 
inadequate  such  as  around  the  right  hand  end  of  the  permanent  islands  in  site  5  (Figure  5.4e). 
Potentially  this  type  of  information  allows  specific  improvements  in  the  model’s 
bathymetry/topography  to  be  made.  The  effect  of  varying  the  friction  values  on  the  spatial 
predictions  varies  from  site  to  site.  For  example  at  sites  1  and  3  there  is  very  little  variation  in 
inundated  area  but  at  sites  2  and  6  the  size  of  the  island(s)  decreases  substantially  as  the  friction 
increases,  i.e.  as  the  water  level  rises.  This  allows  calibration  of  the  model  to  be  carried  out 
with  such  data  by  varying  bed  friction  until  there  is  a  close  spatial  match  between  observed  and 
predicted.  This  variable  behaviour  is  caused  by  the  gradient  of  the  topography  around  the  water 
surface.  If  the  topography  has  a  shallow  gradient  then  changes  in  the  bed  friction,  changing  the 
water  depth,  causes  changes  in  the  inundated  areas.  Where  the  topographic  gradient  is  steeper 
then  the  same  changes  in  water  depth  shall  have  virtually  no  impact  on  inundated  area.  A 
problem  arises  in  this  specific  case,  as  far  as  calibrating  the  model  on  this  data,  because  the 
topography  above  the  water  surface  was  only  estimated  from  US  Department  of  the  Interior 
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maps  and  is  therefore  probably  not  of  the  accuracy  required  to  enable  much  confidence  to  be 
placed  in  the  predictions  of  inundated  area  on  islands.  However  estimated  parameter  values  to 
calibrate  the  inundated  area  at  sites  2  and  6  are  0.025  and  0.02. 
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Figure  5-4a  Comparison  of  Landsat  TM  image  and  4  alternative  model  predictions  with 
different  bed  friction  (Manning’s  ‘n’)  values  at  site  1  for  event  1. 
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n  =  0.03 
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Figure  5-4b  Comparison  of  Landsat  TM  image  and  4  alternative  model  predictions  with 
different  bed  friction  (Manning’s  ‘n’)  values  at  site  2  for  event  1. 
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Figure  5-4c  Comparison  of  Landsat  TM  image  and  4  alternative  model  predictions  with 
different  bed  friction  (Manning’s  ‘n’)  values  at  site  3  for  event  1. 
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n  =  0.01  Model  Predictions  -  0.02 


n  =  0.03  n  =  0.04 


Figure  5-4d  Comparison  of  Landsat  TM  image  and  4  alternative  model  predictions  with 
different  bed  friction  (Manning’s  ‘n’)  values  at  site  4  for  event  1. 
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Figure  5-4e  Comparison  of  Landsat  TM  image  and  4  alternative  model  predictions  with 
different  bed  friction  (Manning’s  ‘n’)  values  at  site  5  for  event  1. 
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Figure  5-4f  Comparison  of  Landsat  TM  image  and  4  alternative  model  predictions  with 
different  bed  friction  (Manning’s  ‘n’)  values  at  site  6  for  event  1. 
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By  geo-referencing  the  model  results  and  satellite  image  an  overlay  of  the  two  can  be  made.  An 
example  of  such  is  shown  in  Figure  5.5  for  a  section  of  this  reach.  This  detailed  comparison 
between  observed  and  predicted  flow  field  boundary  enables  problem  areas  to  be  pinpointed 
and  further  model  development  to  proceed  by  further  data  collection  or  adjustments  to  the 
model  in  these  areas. 


Figure  5-5  An  overlay  of  the  model  prediction  of  flow  field  boundary  onto  the  Landsat 
TM  image  for  a  region  around  Goat  Island  (including  sites  5  and  6)  for  event  1. 


5.1.3  Summary  of  Event  1 

These  first  simxilated  flow  conditions  have  shown  the  both  the  excellent  predictions  of  the  2-D 
hydraulic  model  and  the  utility  of  the  data  sources  to  validate  its  performance.  The  internal 
stage  data  was  matched  at  single  points  but  not  so  closely  at  both  points  simultaneously.  The 
inundation  extent  was  validated  against  the  satellite  images  very  well.  Some  areas  showed  up 
problems  in  the  models  bathymetry/topography.  Given  possible  future  refinements  in  the  model 
set-up  the  internal  stage  data  should  be  able  to  be  matched  simultaneously  at  both  sites  quite 
easily.  Matching  all  the  flow  field  boundaries  shall  be  more  difiScult  given  their  dependence  on 
the  model  bathymetry.  This  is  however  expected  given  the  relative  strength  of  the  two 
validation  types  as  discussed  earlier. 


5.2  Event  2  (25th  June  1994) 

This  section  of  the  1994  flow  record  was  chosen  to  be  another  period  of  steady  state  flow  but 
with  different  flowrates  than  those  in  the  previous  event  and  temporally  close  to  the  previous 
event  to  eliminate  potential  effects  of  changing  bed  morphology.  A  section  of  the  record  around 
Jime  25th  1994  was  identified  as  fitting  these  criteria.  The  inflow  from  Gavins  Point  Dam  and 
the  James  River  were  virtually  constant  at  795m3/s  and  52m3/s  respectively  (and  the  outflow 
stage  at  Maskell  constant  at  342.73m).  The  internal  stage  data  was  available  at  both  sites  for 
this  time  but  unfortunately  no  satellite  image  could  be  obtained. 

For  these  flow  conditions  the  same  series  of  seven  runs  with  varying  friction  were  carried  out  as 
for  the  calibration  period.  Again  the  results  showed  strong  linear  tendencies  and  have  been 
treated  in  the  same  manner.  In  this  case  the  equations  relating  water  surface  error  to  Manning’s 
’n'  are: 

Yankton  Error  =  -0.327  +  32.44  n 
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Gayville  Error  =  -0.981  +  30.00  n 

Alternatively  relating  water  depth  to  bed  roughness  (n)  thus: 
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Yankton  Water  Depth  =  4.77  +  32.44  n 
Gayville  Water  Depth  ==  2.81  +30.00  n 

As  with  the  regression  lines  for  the  previous  case  the  residuals  with  these  are  very  small 
indicating  a  high  degree  of  confidence  in  the  relationships  within  the  range  tested. 

Again  the  values  of  bed  roughness  were  calculated  that  eliminated  the  error  at  the  two  gauging 
stations  in  turn.  These  turned  out  to  be  0.0101  to  eliminate  the  error  at  Yankton  and  0.0327  to 
eliminate  the  error  at  Gayville.  Whilst  the  stage  at  Yankton  is  predicted  correctly  the  prediction 
at  Gayville  is  0.68m  too  low  (21.84%  of  the  expected  water  depth)  and  conversely  whilst  the 
stage  at  Gayville  is  predicted  correctly  the  prediction  at  Yankton  is  0.73m  too  high  (12.52%). 
Again  a  central  value  of  bed  roughness  can  be  found  that  approximately  equalises  the  errors  at 
the  two  stations.  In  this  case  that  value  is  0.0214  and  the  error  at  Yankton  and  Gayville  are 
0.36m  and  -0.34m  (6.59%  and  9.85%)  respectively. 

The  results  from  this  flow  event  are  very  similar  to  those  for  event  1  despite  the  different  flow 
conditions.  This  shows  that  the  model  has  the  potential  to  be  robustly  calibrated  so  that  the  a 
single  parameter  set  could  be  used  for  any  flow  event  if  the  bed  morphology  remains 
unchanged  as  it  did  between  these  first  two  flow  events. 


5.3  Event  3  (26th  June  1984) 

This  event  is  during  a  period  of  very  high  flows  in  all  the  tributaries  of  the  Missouri  River,  such 
as  the  James  River.  The  volume  of  water  being  released  from  Gavins  Point  Dam  is  therefore 
very  much  reduced  to  counteract  this  and  limit  downstream  flooding.  No  out  of  bank  flow 
occurs  along  the  modelled  reach.  A  Landsat  TM  image  of  the  reach  on  26th  June  1984  confirms 
this.  Again  this  is  a  steady  state  run  but  several  problems  exist  with  the  flow  data: 

•  The  James  River  gauge  at  Scotland  stops  recording  as  the  flowrate  tops  500m3/s  just  before 
the  date  of  this  simulation. 

•  No  data  is  available  for  the  downstream  boundary  at  Masked 

•  No  internal  stage  data  is  available  at  Yankton  and  only  daily  values  at  Gayville. 

The  simulations  are  therefore  more  of  an  approximation  than  usual  but  the  comparison  to  the 
satellite  image  makes  them  worthwhile.  Again  a  full  suite  of  seven  simulations  has  been  carried 
out  with  varying  friction  values. 

The  predicted  water  level  at  Gayville  again  showed  a  very  strong  linear  tendency  across  the  full 
range  of  bed  friction  values  used.  The  regression  equations  for  this  can  be  written  thus: 

Gayville  Error  =  -  1. 384  +  29.45  n 
Gayville  Water  Depth  =  2.73  +  29.45  n 

The  residuals  with  these  regression  equations  are  very  small  once  again  showing  that  a  high 
degree  of  confidence  can  be  placed  in  the  predictions.  The  value  of  bed  friction  needed  to 
eliminate  the  error  in  the  prediction  is  however  0.047 1  outside  the  range  used  to  formulate  the 
relationship.  The  model  was  therefore  run  at  this  value  to  confirm  that  this  extrapolated 
prediction  holds.  The  result  from  this  simulation  shows  that  the  error  is  0.025m  at  Gayville, 
indeed  a  very  low  error  but  not  perfect,  hence  the  relationship  might  be  starting  to  break  down 

4l 


Missouri  River  Modelling  Proierr. 


but  it  is  unclear.  This  calibrated  'n'  value  is  far  higher  than  any  other  in  any  of  the  simulations. 
Two  possible  reasons  for  this  are: 

•  Uncertainties  in  flow  data  mean  that  the  model  could  be  simulating  a  non-physical  system 
hence  needing  an  unrealistically  large  bed  friction  to  compensate. 

•  This  event  was  10  years  before  the  bathymetric  data  was  collected  hence  substantial 
changes  in  the  bathymetry  over  this  period  could  have  occurred  again  causing  the 
substantial  re-calibration  to  be  necessary. 

These  two  and  other  factors  probably  combined  to  cause  this  effect.  It  does  however  illustrate 
the  way  that  model  parameters  can  be  used  to  overcome  deficiencies  in  the  model  set  up. 

The  model  predictions  are  again  compared  to  the  satellite  image  at  the  6  sites  of  interest  (Figure 
5.6a-f).  At  the  first  two  sites,  where  the  flow  is  very  low,  the  bathymetric  representation  in  the 
model  is  shown  to  be  inadequate  except  with  very  high  bed  friction.  A  friction  value  of  about 
0.04  or  higher  is  needed  to  calibrate  the  model  at  these  sites.  At  sites  3,4  and  5  the  comparison 
is  quite  close  and  varies  very  little  with  bed  friction.  At  site  6  the  inundation  does  vary 
considerably  with  friction,  the  closest  match  being  at  0.03,  but  the  positioning  of  the  island 
could  be  much  better.  The  model  can  therefore  be  said  to  predict  the  inundation  well 
considering  the  time  elapsed  between  this  event  and  the  collection  of  the  bathymetric  data. 
Some  problems  were  apparent  with  the  very  low  flow  at  the  top  of  the  reach. 

Furthermore  the  conditions  of  the  flow  in  this  simulation,  the  channel  geometry  and  the 
satellite  image  combine  to  allow  an  unrivalled  degree  of  process  validation  to  be  performed. 
Process  validation  is,  as  the  name  suggests,  simply  comparing  the  processes  in  the  field  and  the 
model.  Potentially  it  is  the  most  powerful  form  of  model  validation  as  it  addresses  the 
fundamentals  of  physically  based  modelling  rather  than  simply  the  end  results.  The  three 
factors  allowing  this  form  of  validation  to  be  carried  out  in  this  specific  case  are: 

•  The  flow  in  the  main  stem  is  low  but  in  the  James  River  it  is  very  high,  more  than  likely  out 
of  bank  further  up  the  reach. 

•  The  flood  flows  in  the  James  River  carry  a  large  sediment  load  compared  to  the  clear  water 
released  from  Gavins  Point  Dam  on  the  Missouri  main  stem.  This  facilitates  the 
differentiation  of  the  two  water  sources  as  water  carrying  large  quantities  of  sediment  has  a 
higher  reflectance  than  clear  water  in  the  visible  region  of  the  spectrum  (Lillesand  and 
Kiefer,  1987).  The  Landsat  TM  image  used  for  this  comparison  has  two  wave  bands  in  the 
visible  spectrum  allowing  this  difference  to  be  picked  up. 

•  The  confluence  of  the  two  rivers  occurs  at  the  James  River  Island.  The  flow  from  the  James 
River,  under  normal  conditions,  flows  around  the  north  of  the  island  in  a  narrow  channel 
before  joining  the  main  flow.  Under  flood  conditions  the  possibility  does  exist  for  some  of 
the  flow  to  go  round  the  Island  the  other  way,  joining  the  main  stem  much  sooner. 
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Landsat  TM  image 
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n  =  0.01  n  =  0.02 


Model  Predictions 


n=  0.03  n  =  0.04 


Figure  5.6b  Comparison  of  Landsat  TM  image  and  4  alternative  model  predictions  with 
different  bed  friction  (Manning’s  ‘n’)  values  at  site  2  for  event  3. 
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n  =  0.03 


n  =  0.04 


Figure  5.6c  Comparison  of  Landsat  TM  image  and  4  alternative  model  predictions  with 
different  bed  friction  (Manning’s  ‘n’)  values  at  site  3  for  event  3. 
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Model  Prediction  -  0-02 


n  =  0.03  n  =  0.04 


Figure  5.6d  Comparison  of  Landsat  TM  image  and  4  alternative  model  predictions  with 
different  bed  friction  (Manning’s  ‘n’)  values  at  site  4  for  event  3. 
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n  =  0.01 


n  =  0.02 


Model  Predictions 


Figure  5.6e  Comparison  of  Landsat  TM  image  and  4  alternative  model  predictions  with 
different  bed  friction  (Manning’s  ‘n’)  values  at  site  5  for  event  3. 
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n  =  0.01 


n  =  0.02 


Model  Predictions 


n  =  0.03 


n  =  0.04 


Figure  5.6f  Comparison  of  Landsat  TM  image  and  4  alternative  model  predictions  vt'ith 
different  bed  friction  (Manning’s  ‘n’)  values  at  site  6  for  event  3. 
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Looking  first  at  the  satellite  image  of  the  confluence  (Figure  5.7a)  it  can  clearly  be  seen  that  the 
flow  does  split  around  the  James  River  Island  as  was  hypothesised  as  occurring  in  flood 
conditions.  This  flow  pattern,  around  the  island,  compares  very  well  to  the  model  prediction 
(Figure  5.7b).  The  model  can  therefore  be  said  to  represent  the  flow  processes  occurring  under 
these  novel  conditions. 


This  flow  event  had  substantial  problems  in  its  set-up  but  despite  this  valuable  results  have 
been  achieved.  Firstly  the  internal  stage  data  was  shown  that  it  could  be  matched  even  if  outside 
the  range  of  parameter  values  expected.  Secondly,  flow  field  validation  against  satellite  data 
was  shown  to  be  good  but  some  problems  exist  perhaps  because  of  channel  morphology 
changes  between  when  the  image  was  taken  in  1984  and  when  the  bathymetric  data  was 
collected  in  1995.  Finally  an  example  of  the  most  powerful  form  of  validation,  process 
validation  is  shown.  This  is  a  very  limited  example  of  this  but  does  illustrate  that  the  model 
simulates  the  essential  flow  processes  occurring  in  the  reach. 
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5.4  Event  4  (8th  February  1991) 

Contrary  to  the  previous  event  this  one  is  under  conditions  of  very  low  flows  from  both  the 
Missouri  main  stem  and  virtually  non-existent  flow  in  the  James  River.  This  creates  a  new  set 
of  problems  for  the  model.  A  full  set  of  internal  data  and  a  high  resolution  SPOT  image  taken 
on  8th  February  1 99 1  of  most  of  the  reach  are  available. 

A  brief  inspection  of  the  image  shows  that  the  low  flows  cause  the  channel  to  become  braided, 
the  flow  fdling  only  a  tiny  portion  of  its  usual  channel.  Comparison  of  this  to  the  first  model 
run  (n  =  0.025)  shows  instantly  that  the  model  is  not  highly  resolved  enough  or  containing  great 
enough  bathymetric  data  to  allow  the  simulation  of  the  braided  channel  (Figure  5.8).  Further 
runs  would  therefore  be  irrelevant. 

A  quick  comparison  of  the  predicted  to  observed  stage  values  at  Yankton  and  Gayville  is 
shown  in  Table  5.2.  The  match  is  reasonable  but  could  no  doubt  be  calibrated  to  be  much  closer 
but  again  this  is  pointless  given  the  evidence  showing  the  model  is  failing  to  represent  the 
processes  under  these  flow  conditions. 


Table  5-2  Comparison  of  observed  and  predicted  water  surface  elevations  at  Yankton  and 
Gayville  on  8th  February  1991. 


Yankton 

Gayville 

Observed  (m) 

"  ■  ,  ’  35l;l3  . 

:  348.24, 

Predicted  (m) 

■  351.58  . 

,348.84 

Simulating  this  flow  event  has  highlighted  problems  in  the  model  performance  under  very  low 
flow  conditions.  The  model  bathymetry  is  not  refined  enough  to  simulate  the  braided  channel 
that  forms  under  these  conditions  as  shown  on  the  satellite  image. 
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Figure  5.8  Comparison  of  SPOT  image  to  model  predictions  for  event  4  at  (c)  site  4  and  (b) 
site  6. 
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Figure  5.8  Comparison  of  SPOT  image  to  model  predictions  for  event  4  at  (e)  site  5. 
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5.5  Event  5  (13th  -  14th  June  1994) 

Previous  events  have  all  been  modelled  under  steady  state  conditions  but  this  event  is  the  one 
dynamic  section  of  the  flow  record  that  has  been  simulated.  Hence  simulation  of  this  event  is 
the  first  step  in  analysing  the  performance  of  Missouri  River  application  of  TELEMAC-2D 
with  conditions  varying  with  time.  The  event  covers  a  twelve  hour  period  from  13.00  hrs  on  the 
13th  June  1994  until  01.00  hrs  on  the  14th.  During  this  period  the  inflow  from  the  James  River 
remains  constant  at  45m3/s  but  the  inflow  from  Gavins  Point  Dam  decreases  from  905mVs  to 
ISSvnls  between  14.00  hrs  and  15.00  hrs  on  the  13th.  The  outflow  elevation  at  Masked 
decreases  from  342.9m  to  342.8m  during  the  final  6  hours  of  the  simulation.  Initial  conditions 
for  the  simulation  are  the  steady  state  using  these  starting  values,  which  are  that  exact 
conditions  for  event  1.  No  satellite  image  is  available  for  this  event  but  stage  elevations  at  both 
Yankton  and  Gayville  are  available  hourly  throughout.  Performance  can  therefore  be  analysed 
through  time  but  only  at  two  points  in  space. 

Simulations  of  this  event  were  carried  out  for  Manning's  'n'  values  from  0.01  to  0.04  in 
increments  of  0.005,  the  predicted  and  observed  results  at  Yankton  and  Gayville  are  plotted  in 
Figure  5.9.  Looking  first  at  the  behaviour  of  the  model  over  this  dynamic  event  shows  that  the 
bed  friction  influences  not  only  the  water  surface  elevation,  as  seen  before,  but  also  the  timing 
of  the  changes  at  the  gauge  stations.  This  has  potential  implications  to  the  calibration  of  the 
model.  Now  the  predicted  and  observed  data  are  compared  at  the  two  stations. 

With  the  initial  conditions  of  these  simulations  being  identical  to  those  in  event  1  the  calibrated 
'n'  values  from  that  event  can  be  used  here.  Hence  initially  the  calibrated  value  at  Yankton  is 
0.01 18  and  Gayville  is  0.0306.  Figure  5.9(a  and  b)  confirms  this.  At  the  end  of  the  simulations 
the  stage  values  (observed  and  predicted)  at  Yankton  have  returned  to  a  steady  state.  At 
Gayville,  further  down  the  reach,  they  are  still  changing.  Immediately  noticeable  from  Figure 
5.9a  is  that  at  the  end  of  the  simulation  the  calibrated  'n'  value  at  Yankton  is  different  than  at 
the  start.  The  same  is  likely  true  at  Gayville  if  Figure  5.9b  is  extrapolated.  At  Yankton  the  ideal 
parameter  value  decreases  over  the  simulation  but  at  Gayville  it  seems  to  increase.  This  is 
consistent  with  event  2  where  the  flow  conditions  are  similar  to  those  at  the  end  of  this  event. 
In  effect  this  shows  that  at  Yankton  the  model  predicts  a  larger  change  in  water  surface 
elevation  than  is  observed  and  at  Gayville  the  opposite  occurs  (Table  5.3).  The  nature  of  the 
change  of  the  stage  values  is  slightly  awkward  to  analyse  given  the  relative  low  density  of  the 
observed  data.  Figure  5.9a  shows  that  the  predicted  decline  in  water  surface  elevation  at 
Yankton  is  earlier  and  quicker  but  not  as  great  as  the  observed.  At  Gayville  (Figure  5.9b)  the 
predicted  reduction  is  later,  slower  and  larger  than  the  observed. 


Table  5-3  Observed  and  predicted  changes  in  water  surface  elevation  over  dynamic  event 
5. 


Yankton 

Gayville 

Observed  Change  (m) 

Predicted  Change  (m) 

-  -  0.140  1 

These  simulations  on  the  Missouri  River  application  of  TELEMAC-2D  show  that  the  model 
does  predict  the  change  in  water  surface  elevation  at  two  points  fairly  accurately  under 
dynamic  stresses,  the  timing  and  magnitude  of  the  changes  being  broadly  comparable  with  the 
observed  given  the  overall  uncertainties  involved.  Two  problems  do  however  stand  out; 

•  The  continuing  different  behaviour  at  the  two  gauge  stations,  this  time  extending  over  the 
full  dynamic  event. 
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•  Only  having  a  single  parameter  to  vary  means  that  the  water  surface  elevation  and  timing 
cannot  be  varied  separately  but  change  together  with  friction. 

Overall  however  this  initial  test  shows  that  this  model  is  capable  of  simulating  dynamic  flow 
events  with  a  high  degree  of  accuracy  within  the  model  domain.  Further  tests,  preferably  with 
spatially  distributed  comparison  data,  would  be  necessary  to  make  a  more  comprehensive  case 
for  the  models  application  in  such  events. 
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Figure  5-9  Observed  and  predicted,  for  several  bed  friction  values,  water  surface 
elevations  at  (a)  Yankton  and  (b)  Gayville  over  the  duration  of  dynamic  event  5. 
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5.6  Summary  of  the  comparison  to  observed  data 

The  model  has  been  tested  over  widely  varying  flow  conditions  and  compared  to  data  through 
time  and  space.  In  general  the  model  has  been  shown  to  be  a  good  simulator  of  the  river 
hydraulics  under  most  conditions.  The  only  serious  problems  occurred  under  exceptionally  low 
flow  conditions.  The  available  data  allowed  detailed  comparisons  of  observed  and  predicted 
giving  an  invaluable  portrayal  of  model  performance.  The  strength  of  the  model  validation  has 
also  been  assessed.  The  weakest  form  of  validation  here  is  the  internal  stage  data.  This  has  been 
matched  to  about  0.3m  at  point  sites  on  a  first  approximation  of  the  model  set-up.  Further 
changes  to  the  model  parameters  could  reduce  these  errors  significantly.  The  next  step  up  in 
validation  is  the  2-D  flow  field  boundary  validation  against  the  satellite  imagery.  Again  this 
was  shown  to  be  good  for  a  first  attempt.  Some  areas  in  the  model  bathymetry  were  shown  to 
be  inadequate  to  meet  this  level  of  validation.  Further  data  collection  and  model  refinement 
would  be  necessary  to  meet  this  higher  level  of  model  validation.  A  simple  example  of  process 
validation  was  also  illustrated.  The  model  was  shown  to  represent  the  essential  flow  processes 
in  this  particular  limited  example. 

The  first  three  objectives  of  this  study  have  therefore  now  been  achieved.  A  high-space  time 
distributed  model  of  river  hydraulics  has  been  produced  and  shown  to  be  a  good  representation 
of  the  fluvial  system.  The  validation  data  supplied  by  the  Corps  of  Engineers  has  also  been 
shown  to  be  invaluable  to  the  assessment  of  model  behaviour  and  its  effectiveness  in  achieving 
an  improved  partial  validation  of  the  model  estimated. 
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6.  Mesh  Resolution  Comparison 

The  two  meshes  and  their  relative  sensitivities  have  been  discussed  previously.  This  section 
looks  first  at  some  aspects  of  the  differing  bathymetric  representation  between  the  two  meshes 
and  then  the  predictions  of  the  high  resolution  mesh  are  compared  to  the  observed  data  for  the 
steady  state  event  1.  Again  simulations  have  been  carried  out  for  varying  bed  fiiction  values 
fi-om  0.01  to  0.04,  incrementing  by  0.005. 

The  higher  resolution  mesh  allows  a  more  detailed  representation  of  the  bathymetry  of  the  river. 
However  whether  this  representation  is  more  accurate  is  unclear  as  the  bathymetric  data  is  more 
widely  dispersed  than  even  the  low  resolution  mesh.  Figure  6. 1  shows  an  example  bathymetric 
representation  of  the  two  meshes  at  one  areas  of  the  reach.  Although  the  bathymetries  are 
broadly  similar,  as  would  be  expected,  small  differences  do  exist  around  potentially  sensitive 
areas  such  as  sand  banks  or  channel  boundaries.  Whether  they  are  improved  in  the  high 
resolution  mesh  shall  be  explored  later  when  predictions  are  compared  to  the  satellite  image. 
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Figure  6-1  Comparison  of  a  section  of  the  model  bathymetries  in  the  (a)  low  and  (b)  high 
resolution  meshes. 


Missouri  River  Modelling  Project 


6.1  Comparison  to  observed  data 

Firstly,  however,  the  performance  of  the  higher  resolution  model  at  the  two  internal  gauge 
stations  is  examined.  Regression  equations  have  again  been  formulated  as  the  models  response 
is  linear  at  both  stations  and  can  be  written  thus: 

Yankton  Water  Depth  =  4.39  +  32.63  n 
Gayville  Water  Depth  =  2.24  +  33.00  n 

The  results  from  applying  these  equations  and  the  comparison  with  the  lower  resolution  mesh 
can  be  seen  in  Table  6.1.  From  this  table  of  results  it  can  be  seen  that  the  bed  friction  values  are 
consistently  larger  with  the  higher  resolution  mesh  and  the  errors  are  also  greater  in  ever  case. 
Using  the  more  highly  specified  model  cannot  therefore  be  justified  on  this  evidence. 


Table  6-1  Comparison  of  the  calibrated  results  from  the  two  meshes. 


Mesh  1 

Mesh  2 

'n'  value  for  zero  error  at  Yankton 
Corresponding  error  at  Gayville  (m) 

0.0124 

— O.T^:(28.69%)  '  - 

'n'  value  for  zero  error  at  Gayville 
Corresponding  error  at  Yankton  (m) 

*<!  ■  0,0306  •  ' 

6.65  (10.96%)  .  1 ' 

'  ^;0354  - 

bi7lb3.5#o)' 

Intermediate  value  of 'n' 
Corresponding  error  at  Yankton  (m) 
Corresponding  error  at  Gayville  (m) 

.  0.0212  '  : 
'  0.33  (5.89%) . 
-0.30  (8.37%) 

The  same  six  sites  have  been  chosen  for  the  spatial  comparison  of  model  results  from  the  high 
resolution  mesh  to  the  Landsat  TM  image  as  were  used  earlier  with  the  lower  resolution  mesh. 
Again  the  images  are  compared  side  by  side  to  four  model  predictions  with  varying  bed  friction 
(Figure  6.2).  It  is  immediately  evident  that  the  resolution  of  the  predictions  is  higher,  than  for 
the  low  resolution  mesh  (Figure  5.4)  but  whether  this  offers  an  improvement  in  prediction  is 
however  not  so  clear.  Closer  inspection  of  the  figures  suggest  that  in  some  areas  the  new 
predictions  do  have  a  closer  spatial  match,  such  as  at  sites  3  and  4,  but  this  is  by  no  means 
conclusive. 

Compared  to  the  low  resolution  mesh,  figure  6.2  shows  that  to  attain  the  same  amount  of 
inundation  the  friction  values  must  be  higher  with  the  high  resolution  mesh.  This  is  evident  at 
sites  1,2,4  and  6  and  concurs  with  the  results  for  the  stage  values.  Therefore  it  can  be  said  that 
to  achieve  the  same  water  surface  elevation  requires  a  higher  friction  value  in  the  high 
resolution  mesh.  Why  this  should  be  so  is  unclear  but  could  be  a  result  of  differing  topographic 
or  parameter  representations. 
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Figure  6-2a  Comparison  of  Landsat  TM  image  to  4  alternative  model  predictions  with 
different  values  of  bed  friction  (Manning’s  ‘n’)  on  the  high  resolution  mesh  (mesh  2)  for 
event  1  at  site  1. 
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Figure  6-2b  Comparison  of  Landsat  TM  image  to  4  alternative  model  predictions  with 
different  values  of  bed  friction  (Manning’s  ‘n’)  on  the  high  resolution  mesh  (mesh  2)  for 
event  1  at  site  2. 
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Figure  6-2c  Comparison  of  Landsat  TM  image  to  4  alternative  model  predictions  with 
different  values  of  bed  friction  (Manning’s  ‘n’)  on  the  high  resolution  mesh  (mesh  2)  for 
event  1  at  site  3. 
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Landsat  TM  image 


1000m 


n  =  0.03  n  =  0.04 


Figure  6-2d  Comparison  of  Landsat  TM  image  to  4  alternative  model  predictions  with 
different  values  of  bed  friction  (Manning’s  ‘n’)  on  the  high  resolution  mesh  (mesh  2)  for 
event  1  at  site  4. 
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Figure  6-2e  Comparison  of  Landsat  TM  image  to  4  alternative  model  predictions  with 
different  values  of  bed  friction  (Manning’s  ‘n’)  on  the  high  resolution  mesh  (mesh  2)  for 
event  1  at  site  5. 
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Figure  6-2f  Comparison  of  Landsat  TM  image  to  4  alternative  model  predictions  with 
different  values  of  bed  friction  (Manning’s  ‘n’)  on  the  high  resolution  mesh  (mesh  2)  for 
event  1  at  site  6. 
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6.2  Summary 

This  brief  look  at  the  impact  of  mesh  resolution  on  the  model  results  has  shown  this  to  be  an 
important  factor.  The  bed  friction  values  must  be  higher  to  achieve  the  same  results  than  with 
the  lower  resolution  mesh  and  the  bathymetry  is  of  a  higher  resolution,  if  not  more  accurate. 
Mesh  resolution  is  a  continual  problem  with  no  obvious  solution  except  perhaps  a  compromise 
between  increased  representation  with  a  high  resolution  mesh  and  reduced  computational 
expense  with  a  lower  resolution  mesh.  Where  this  compromise  lies  is  difficult  to  judge  unless 
increased  validation  data  is  available  to  judge  model  performance  in  even  more  detail  to  justify 
increasing  the  representation  in  terms  of  model  results  improvements.  From  this  study  this 
would  not  seem  to  be  the  case  given  the  limitations  imposed  by  the  data,  particularly 
bathymetry.  The  following  section  continues  to  look  at  the  bathymetric  data  in  terms  of  its 
impact  on  model  predictions. 
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7.  Bathymetric  Data  Requirements 

The  high  cost  of  collecting  bathymetric  and  topographic  data  must  be  justified  by  attaining  high 
quality  model  results  that  could  not  be  accomplished  with  a  more  basic  representation.  The 
impact  of  the  density  of  cross  section  on  the  models  bathymetry  and  model  predictions  are 
investigated  in  this  section. 

Four  realisations  of  the  model  bathymetry  were  created  using: 

•  all  the  available  cross  sections  (150  m  spacing), 

•  alternate  cross  sections  (300  m  spacing), 

•  one  in  four  cross  sections  (600  m  spacing), 

•  one  in  eight  cross  sections  (1200  m  spacing). 

From  each  of  these  data  sets  the  model  bathymetry  was  created  at  the  node  points  on  the  low 
resolution  mesh  using  the  directional  interpolator  described  in  section  3.2.  The  model  was  then 
run  once  with  each  of  these  representations  for  the  conditions  of  event  1  with  a  bed  fiiction 
value  of  0.025. 

The  impact  of  these  changes  on  the  model  bathymetry  was  profound.  Figure  7.1  shows  the 
impacts  in  one  example  area  (site  6).  As  can  be  clearly  seen  the  detail  of  the  bathymetric 
representation  in  the  model  degrades  rapidly  as  the  number  of  cross  sections  decrease. 
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Figure  7-1  Impact  of  degrading  the  number  of  cross  sections  used  on  the  model 
bathymetry.  The  spacing  degrades  from  150  m  to  1200  m  across  the  range  illustrated. 
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7.1  Comparison  to  observed  data 

The  influence  of  changing  the  level  of  bathymetric  data  on  the  water  surface  elevation 
predictions  at  Yankton  and  Gayville  is  shown  in  Table  7. 1 .  From  this  table  it  can  be  seen  that  at 
Yankton  the  level  of  data  provision  has  minimal  effects  of  the  predictions,  varying  them  by 
only  O.lm  compared  to  an  actual  error  of  approaching  0.5m.  At  Gayville,  however,  decreasing 
the  bathymetric  produces  a  continuous  increase  in  water  surface  elevation.  The  magnitude  of 
the  change  over  the  runs  is  0.72m,  spanning  the  observed  value.  The  half  data  set  run  produces 
a  prediction  of  water  surface  elevation  at  Gayville  within  0.015m  of  the  observed  value. 


Table  7-1  Water  surface  elevations  at  Yankton  and  Gayville  for  simulations  with  the 
various  levels  of  bathymetric  data  supplied  to  the  model. 


Location 

Full 

Half 

Quarter 

Eighth 

Observed 

Yankton 

352.667 

::;352.^TS;f|||:, 

352.196 

Gayville 

349.854 

349.990 

350.398 

350.032 

From  these  results  there  appears  to  be  little  benefit  in  using  the  highest  density  data  set  as  the 
results  are  equally  good  using  half  the  data.  The  quarter  and  eighth  runs  do  however  show 
decay  in  the  quality  of  the  results  at  Gayville.  The  more  important  results  are  however  the 
spatially  distributed  ones  that  are  now  investigated. 

Once  more  the  side  by  side  comparisons  have  been  made  for  model  predictions,  this  time  with 
the  4  levels  of  bathymetric  data,  against  the  Landsat  TM  image  at  the  6  sites  of  interest  (Figure 
7.2).  In  each  case  the  full  bathymetry  appears  to  give  a  favourable  comparison  against  the 
observed  data.  Reducing  the  bathymetric  data  has  a  major  influence  on  the  spatial  predictions. 
In  some  cases  this  manifests  itself  even  when  half  the  data  is  used  such  as  at  sites  1,3  and  6, 
where  the  spatial  match  degrades  rapidly.  As  the  data  reduction  continues  huge  errors  are 
introduced  at  some  of  the  sites,  e.g.  site  1,  but  not  so  in  other,  e.g.  site  2. 


7.2  Summary 

The  results  from  this  section  have  shown  the  bathymetric  data  to  have  a  profound  influence  on 
the  spatial  predictions  of  the  2-D  hydraulic  model  and  a  lesser  influence  on  the  water  surface 
elevation.  This  fits  well  with  the  rigour  of  validation  discussed  in  section  5.  There  it  was 
proposed  that  matching  internal  stage  values  was  a  weaker  form  of  validation  than  matching  2- 
D  inundation.  These  results  show  that  a  more  poorly  specified  model,  i.e.  that  with  less 
bathymetric  data,  can  match  the  internal  stage  values  but  not  the  inundation  extent.  Flence 
inundation  extent  validation  is  a  more  rigorous  test  of  model  performance. 

From  these  results  it  would  seem,  as  expected,  that  the  highest  density  of  bathymetric  data 
produces  the  best  results.  More  detailed  simulations  looking  at  the  relationship  between  mesh 
resolution,  bathymetric  data  density  and  distribution,  and  interpolation  method  must  be  carried 
out  before  more  detailed  guidelines  on  the  optimum  level  of  bathymetric  data  can  be 
considered. 
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Figure  7-2a  Comparison  of  the  Landsat  TM  image  to  model  predictions  with  the  4  levels  of 
model  bathymetry  used  for  event  1  at  site  1.  The  bathymetries  were  created  by  varying  the 
spacing  of  cross  sections  used  from  150  m  to  1200  m. 
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Figure  7-2b  Comparison  of  the  Landsat  TM  image  to  model  predictions  with  the  4  levels  of 
model  bathymetry  used  for  event  1  at  site  2.  The  bathymetries  were  created  by  varying  the 
spacing  of  cross  sections  used  from  150  m  to  1200  m. 
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Figure  7-2c  Comparison  of  the  Landsat  TM  image  to  model  predictions  with  the  4  levels  of 
model  bathymetry  used  for  event  1  at  site  3.  The  bathymetries  were  created  by  varying  the 
spacing  of  cross  sections  used  from  150  m  to  1200  m. 
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Figure  7-2d  Comparison  of  the  Landsat  TM  image  to  model  predictions  with  the  4  levels  of 
model  bathymetry  used  for  event  1  at  site  4.  The  bathymetries  were  created  by  varying  the 
spacing  of  cross  sections  used  from  150  m  to  1200  m. 
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Figure  7-2e  Comparison  of  the  Landsat  TM  image  to  model  predictions  with  the  4  levels  of 
model  bathymetry  used  for  event  1  at  site  5.  The  bathymetries  were  created  by  varying  the 
spacing  of  cross  sections  used  from  150  m  to  1200  m. 
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Figure  7-2f  Comparison  of  the  Landsat  TM  image  to  model  predictions  with  the  4  levels  of 
model  bathymetry  used  for  event  1  at  site  6.  The  bathymetries  were  created  by  varying  the 
spacing  of  cross  sections  used  from  150  m  to  1200  m. 


Missouri  River  Modelling  Prnierr 


8.  Concluding  Remarks 

This  modelling  study  of  a  55  km  reach  of  the  Missouri  River  has  illustrated  the  capability  of 
high  resolution,  physically  based  fluvial  hydraulic  models  for  accurate  flow  representation.  A 
hydraulic  model  was  created  for  the  reach  of  interest,  between  Gavins  Point  Dam  and  Maskell, 
using  the  TELEMAC-2D  modelling  system.  A  sensitivity  analysis  of  the  model  showed  bed 
friction  to  be  the  dominant  parameter  governing  model  performance.  Extensive  validation  of 
the  model  predictions  was  carried  out  against  both  internal  stage  data  and  satellite  imagery  of 
flow  field  boundaries.  Finally  the  impact  of  mesh  resolution  and  bathymetric  data  provision 
was  briefly  examined. 

The  key  findings  of  this  study  can  be  summarised  thus; 

•  A  state  of  the  art  two-dimensional  hydraulic  model  can  simulate  fluvial  hydraulics 
accurately  on  this  scale  in  both  steady  state  and  dynamic  conditions. 

•  Model  predictions  can  be  matched  to  stage  data  internal  to  the  model  domain  with  a  high 
degree  of  accuracy  for  initial  simulations  (Figure  5.2).  Improvements  to  the  match  could  no 
doubt  follow  with  further  model  developments. 

•  Spatial  distribution  of  flow  field  boundaries  has  been  shown  to  be  accurate  in  most  cases  by 
comparison  to  satellite  imagery  (Landsat  TM  and  SPOT)  (Figure  5.4). 

•  Partial  model  validation  has  been  graded  by  its  rigour.  Matching  2-D  inundation  data  is 
considered  more  powerful  than  matching  internal  stage  data,  a  hypothesis  corroborated  by 
model  predictions  with  varying  specifications  of  bathymetry  (Figure  7.2  and  Table  7.1). 

•  A  degree  of  process  validation  was  achieved  against  satellite  imagery  in  a  specific  case 
(Figure  5.7). 

•  The  relationship  between  mesh  resolution  and  bathymetric  data  is  vital  to  model 
specification  but  requires  further  study  to  quantify. 


There  remain  however,  many  areas  where  further  work  could  improve  the  model  predictions. 

•  The  fundamental  complexity  and  connection  of  mesh  and  bathymetry  is  a  limiting  factor 
whilst  a)  only  limited  bathymetric  data  is  available  and  b)  computational  expense  must  be 
considered  when  making  meshes.  Until  meshes  can  be  much  more  finely  resolved  and  the 
bathymetric  data  available  to  apply  to  them  at  a  commensurate  scale  then  this  problem  of 
bathymetric  representation  and  its  impact  on  model  predictions  shall  remain.  Further  data 
collection  in  specific  problem  areas  could  mitigate  the  problem  in  this  specific  case. 

•  The  above  problem  is  manifested  predominantly  under  low  flow  conditions  (Figure  5.8) 
where  the  bed  representation  is  exposed  more  than  usual. 

•  The  model  set-up  used  here  in  terms  of  parameters  was  extremely  simple  and  as  such 
model  results  were  not  as  good  as  they  could  be  following  further  modifications.  With 
modification  of  the  bed  friction  parameter,  distributing  it  over  time  and/or  space  the  errors 
particularly  in  internal  stage  could  be  reduced  dramatically.  The  2-D  flow  field  boundaries 
are  however  more  dependant  on  the  bathymetry  than  the  parameter  distribution. 

•  With  the  knowledge  gained  from  the  simulations  carried  out  already  the  model  structure 
could  be  modified  to  include  more  permanent  islands  where  they  are  evident  and  vary  the 
mesh  resolution  over  the  domain  where  the  gradient  of  predicted  variables  are  greatest. 

•  Now  that  model  validation  has  been  somewhat  structured  the  possibility  exists  for 
deliberate  data  collection  in  order  to  validate  specific  outputs  of  the  model  to  produce  as 
strong  and  complete  a  validation  as  possible.  For  example  more  2-D  process  data  could  be 
collected  in  order  to  validate  the  most  fundamental  aspect  of  a  2-D  hydraulic  model. 
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The  ever  progressing  capability  of  highly  resolved  2-D  fluvial  hydraulic  models,  including 
flows  over  dry  areas,  opens  up  a  wide  range  of  potential  applications  that  include: 

•  Predicting  the  impact  of  anthropogenic  activity 

•  Flood  hazard  assessment 

•  Geomorphological  and  hydrological  process  studies 

•  River  management  schemes 

This  report  has  shown  that  the  capability  now  exists  for  high  resolution  fluvial  hydraulic 
modelling  to  be  applied  over  long  (>50  km)  reaches  with  a  high  degree  of  accuracy.  The  use  of 
wetting  and  drying  algorithms  allows  the  simulation  of  dynamic  flows  over  initially  dry  areas 
and  their  retreat.  This  capability  opens  up  a  wide  range  of  potential  applications  especially  for 
the  “what  if’  type  modelling  vital  for  controlled  management  and  development  of  river 
systems. 
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